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Abstract 

One  of  the  primary  problems  with  the  application  of  Space-Time  Adaptive  Processing  (STAP) 
techniques  is  secondary  data  support  for  the  interference  plus  noise  covariance  matrix  estimate. 
This  estimate  is  obtained  using  a  Maximum  Likelihood  Estimator  (MLE).  As  always  in  the  ceise 
of  a  statistical  estimator,  the  larger  the  number  of  points,  or  secondary  data  support,  used  in  the 
calculation  the  more  accurate  the  estimate. 

Reed  [15]  has  shown  the  required  secondary  data  support  in  this  estimate  to  achieve  perfor¬ 
mance  within  3  dB  of  the  optimal  signal-to-interference-plus-noise  ratio  (SINR)  is  approximately 
two  times  the  Degrees  Of  Freedom  (DOF).  Hence,  reducing  the  DOF  used  in  the  adaptive  algo¬ 
rithm  also  reduces  required  secondary  data  support.  Reed  developed  this  rule  for  Sample  Matrix 
Inversion  (SMI)  STAP  techniques. 

A  concern  arises  when  applying  this  rule  to  a  class  of  newer  STAP  methods  designed  to 
reduce  the  DOF.  An  example  of  one  of  these  newer  methods  is  the  Cross  Spectral  Metric  (CSM) 
algorithm  [6].  Since  the  CSM  algorithm  does  not  fall  under  the  SMI  classification,  application  of 
Reed’s  rule  for  required  secondary  data  support  is  brought  into  question.  In  an  effort  to  answer  this 
question,  extensive  Monte  Carlo  experiments  were  performed.  These  experiments  have  resulted  in 
a  number  of  interesting  conclusions. 

When  using  full  dimensionality  (all  DOF),  the  CSM  algorithm  performs  as  predicted  by 
Reed’s  rule.  However,  when  reducing  the  DOF  the  required  secondary  data  support  does  not 
follow  Reed’s  prediction.  Optimal  SINR  performance  is  obtained  only  when  the  DOF  are  equal  to 
the  dimension  of  the  interference  subspace.  Furthermore,  there  is  a  threshold  below  the  dimension 
of  the  interference  subspace  where  reducing  the  DOF  does  not  gain  any  reduction  in  secondary 
data  support. 

A  second  goal  of  the  thesis  determines  the  impact  of  non-homogeneities  within  the  secondary 
data  on  the  CSM  algorithm.  The  Generalized  Inner  Product  (GIP)  detection  scheme  is  then  used 
to  excise  these  non-homogeneities  from  the  secondary  data.  The  CSM  algorithm  was  found  to 
be  susceptible  to  non-homogeneities.  The  use  of  the  GIP  successfully  negated  the  impact  on  this 
algorithm. 
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Secondary  Data  Support  and  Non-Homogeneities 


in 

Space-Time  Adaptive  Processing 


1.  Introduction 

1.1  Background 

Clutter  and  jamming  are  serious  problems  for  airborne  radar  engineers.  Historically,  the  elim¬ 
ination  of  clutter  has  been  accomplished  through  very  low  antenna  sidelobes,  an  expensive  and  dif¬ 
ficult  proposition  in  the  case  of  airborne  radars.  Furthermore,  Electronic  Counter-Countermeasure 
(ECCM)  techniques  currently  implemented  in  radars  are  reactionary,  meaning  they  are  relegated 
to  respond  defensively  instead  of  offensively.  Space-Time  Adaptive  Processing  (STAP)  promises 
an  inexpensive  alternative  capable  of  filtering  clutter  and  allowing  an  offensive  approach  to  jam¬ 
ming  by  removing  it  before  the  radar  processes  the  return.  STAP  refers  to  adaptive  filtering  of 
antenna  array  inputs  in  both  the  spatial  and  temporal  domains  to  suppress  the  effects  of  clutter 
and  jamming.  Current  research  shows  it  is  an  effective  approach  that  significantly  improves  the 
performance  of  an  airborne  radar. 

STAP  is  a  digital  filtering  technique  applicable  to  any  antenna  array  and  holds  the  potential 
to  significantly  improve  radar  performance  at  minimal  cost.  The  mission  of  Rome  Laboratory  (the 
sponsor  of  this  research  project)  is  to  provide  cost  effective  solutions  to  the  technological  problems 
of  the  Air  Force.  STAP  fits  perfectly  into  their  mission.  This  research  also  supports  the  broader 
Air  Force  goal  to  bring  effective  weapon  platforms  to  the  battlefield  while  minimizing  cost. 
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1,2  Objective 


Although  there  are  a  variety  of  STAP  techniques  available,  each  with  their  own  strengths 
and  weaknesses,  all  of  these  techniques  rely  on  the  statistics  of  the  interference  and  clutter,  or  the 
interference  covariance  matrix  R.  The  widely  varying  clutter  and  interference  characteristics  of  an 
airborne  environment  effectively  prohibits  prediction  of  these  statistics.  The  interference  covariance 
structure  must  be  estimated  from  secondary  data  available  within  the  radar’s  Coherent  Processing 
Interval  (CPI). 

This  thesis  focuses  on  a  method  of  STAP  just  recently  introduced.  The  Cross  Spectral  Metric 
(CSM)  algorithm  [6]  results  in  an  adaptive  choice  of  the  Degrees  of  Freedom  (DOF)  that  result 
in  optimal  SINK.  The  cross-spectral  metric  approach  was  developed  by  Goldstein  and  Reed  [8]  as 
an  alternative  to  the  method  of  principle  components.  Chapter  IV  discusses  these  approaches  in 
more  detail.  Because  of  the  methods  used  to  determine  the  optimum  DOF,  there  is  some  question 
as  to  the  secondary  data  support  required  to  effectively  estimate  the  covariance  structure  of  the 
environment.  This  research  effort  attempts  to  answer  this  question. 

Once  the  required  secondary  data  support  has  been  determined  for  the  CSM  algorithm,  the 
thesis  progresses  to  another  problem  prevalant  in  all  STAP  algorithms.  If  a  target  or  some  other 
non-homogeneity  is  present  within  the  secondary  data,  the  covariance  matrix  estimate  is  distorted. 
Depending  on  the  severity  of  the  non-homogeneity  in  the  secondary  data,  this  estimate  can  distort 
the  antenna  pattern  produced  by  the  STAP  algorithm.  This  research  project  illustrates  the  impact 
of  secondary  data  non-homogeneities  on  the  CSM  algorithm. 

In  an  effort  to  alleviate  the  effect  of  these  non-homogeneities,  Melvin  [14]  has  proposed  a 
non-homogeneity  detection  scheme.  Detection  is  performed  through  the  use  of  a  Generalized  Inner 
Product  (GIP).  He  has  shown  the  improvements  gained  through  its  use  for  the  reduced  dimension 
Factored  Time-Space  (FTS)  approach  [14].  As  a  final  objective,  this  research  applies  the  GIP 
detection  scheme  to  the  CSM  algorithm.  Matlab®  is  used  to  explore  these  alternatives. 
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1.3  Material  and  Equipment 


All  simulations  were  developed  and  performed  using  Matlab®  versions  4.2c  and  5.0,  a  high 
performance  numeric  computation  and  visualization  software  optimized  for  matrix  operations.  The 
software  was  operated  on  Sun  Workstations  provided  by  the  Air  Force  Institute  of  Technology. 

1.4  Support 

Capt  William  L.  Melvin,  Rome  Laboratory,  Griffis  Industrial  Park,  NY  served  as  the  point  of 
contact  for  this  thesis.  Many  thanks  go  to  Rome  Laboratory  for  allowing  the  use  of  Multi-Channel 
Airborne  Radar  Management  (MCARM)  data.  This  data  represents  real  world  collections  from 
an  antenna  array  using  the  MCARM  airborne  radar  platform.  Readers  interested  in  more  details 
should  reference  the  final  technical  report  [18]. 

1.5  Organization 

This  document  serves  to  evaluate  required  secondary  data  support  for  the  CSM  algorithm  and 
the  impact  of  non-homogeneities  among  secondary  data  used  for  interference  plus  noise  covariance 
matrix  estimation.  In  order  to  present  a  coherent  approach,  the  fundamentals  of  beamforming  and 
STAP  are  first  explored.  Chapter  II  delves  into  the  signal  processing  viewpoint  of  antenna  array 
beamforming.  The  basic  concepts  used  in  the  STAP  development  are  derived  and  discussed  there. 

Chapter  III  develops  the  basics  of  STAP  in  parallel  with  Ward  [21].  Both  fully  adaptive 
and  partially  adaptive  algorithms  are  introduced.  Because  of  the  computational  burden  associated 
with  fully  adaptive  STAP,  this  algorithm  is  beyond  current  computing  capabilities  for  real  time 
operation  in  an  airborne  platform.  Although  there  are  many  reduced  dimension  STAP  algorithms 
available,  this  chapter  only  introduces  the  Factored  Time-Space  (FTS)  method.  The  FTS  method 
is  simple  to  implement  and  effectively  introduces  the  advantages  of  reduced  dimensionality  with 
respect  to  the  computational  load  of  STAP. 
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Chapter  IV  is  devoted  to  developing  an  understanding  of  Goldstein’s  [8]  Cross  Spectral  Metric 
(CSM)  algorithm.  The  CSM  is  introduced  as  a  means  of  selecting  an  optimum  dimension  reduction 
matrix  used  to  reduce  the  DOF  and  indirectly  the  secondary  data  support. 

The  final  set  of  background  information  necessary  is  a  complete  discussion  of  the  Generalized 
Inner  Product  (GIP)  non-homogeneity  detector  [14].  Chapter  V  not  only  introduces  this  detection 
scheme  but  also  presents  a  comprehensive  example  of  its  application  to  the  FTS  algorithm  using 
MCARM  data. 

The  analysis  of  the  CSM  algorithm  is  contained  in  Chapter  VI.  There  are  two  major  questions 
answered  in  this  thesis.  The  first  is  in  regard  to  sample  support  required  for  the  CSM  model.  Monte 
Carlo  simulation  is  employed  to  answer  this  question.  The  second  question  of  interest  involves 
application  of  the  GIP  non-homogeneity  detector  to  the  CSM  model.  The  goal  is  to  determine  if 
performance  is  improved. 
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IL  General  Beamforming  Concepts 


2.1  Introduction 

People  are  bombarded  with  plane  waves  every  day,  but  the  average  person  doesn’t  think  of 
them  in  this  manner.  The  greeting  from  a  friend  in  the  distance  is  carried  on  a  plane  wave,  an 
acoustic  wave  in  this  simple  example.  The  visual  image  of  that  friend  is  carried  on  plane  waves  from 
a  relatively  small  portion  of  the  electromagnetic  spectrum.  Extending  these  interpretations,  radar 
and  sonar  become  simple  extensions  of  the  basic  human  sensory  organs.  They  give  us  the  ability 
to  make  use  of  portions  of  the  electromagnetic  spectrum  beyond  our  own  physiological  limitations. 

This  first  chapter  introduces  a  mathematical  model  for  these  waves  from  basic  electromagnetic 
principles.  It  then  highlights  the  functional  similarities  between  the  expressions  developed  using 
this  model  and  an  expression  for  the  multidimensional  Fourier  transform.  This  similarity  shows 
the  applicability  of  multidimensional  signal  processing  in  the  analysis  of  these  signals,  referred  to 
as  space-time  signals.  Beamforming  and  Space-Time  Adaptive  Processing  (STAP)  are  found  to  be 
ideal  methods  to  ’’analyze”  these  signals  and  extract  the  information  of  interest. 

This  chapter  follows  closely  the  development  by  Dudgeon  and  Mersereau  [5,  pages  289-303]. 
The  chapter  lays  the  groundwork  for  the  later  concepts  of  STAP. 

2.2  Foundations 

The  mathematical  representation  for  a  propagating  plane  wave  from  electromagnetic  funda¬ 
mentals  is  in  terms  of  a  complex  exponential, 

s{x,t)  = 


gi(a;o<-ko-x) 


(2.1) 
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where  the  variable  Uo  represents  the  temporal  frequency  of  the  wave  and  the  vector  designates 
the  direction  of  the  propagating  wave  and  its  spatial  frequency.  The  dot  product  in  the  expression 
represents  the  inner  product  between  the  two  vectors  and  can  also  be  written  as 

kjx.  (2.2) 


where  the  superscript  ^  represents  the  transpose  operation.  The  spatial  frequency  vector  can  be 
described  by 


ko 


(2.3) 


where  is  the  nominal  wavelength  of  the  propagating  wave  (equal  to  where  c  is  the  velocity 
of  propagation  and  fo  is  the  transmitted  frequency)  and  k  is  a  unit  vector  in  the  direction  of  the 
propagating  wave.  The  form  of  Eqn.  (2.1)  results  in  the  definition  of  a  new  parameter  a, 


^  /  sec  \ 
Uo  \ meter/  ’ 


(2.4) 


where  the  units  of  this  new  quantity  are  the  inverse  of  velocity.  Hence,  it  is  justifiably  called 
the  slowness  vector.  The  ©  subscript  depicts  the  slowness  vector  is  in  the  direction  of  interest,  or 
reference  direction.  Generically,  the  slowness  vector  is  written  without  the  subscript.  The  plane 
wave  expression  describing  a  signal  propagating  in  the  direction  of  interest  written  in  terms  of  the 
slowness  vector  then  becomes 


s{x,t)  =  e^‘^o(t-<xo.x)^  (2.5) 

The  slowness  vector  not  only  describes  the  speed  of  the  wave,  1/|q:o|,  but  also  it’s  direction. 
The  plane  wave  expression  can  be  written  generically  as  a  function  of  the  speed  of  the  wave  and 
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direction  by  dropping  the  subscript  on  a  and  w, 


(2.6) 


Examine  the  standard  Fourier  transform  definition  as  shown  by 


/oo 

•OO 

J  —CO 


(2.7) 


The  expression  illustrates  a  conflict  with  the  sign  in  the  exponential  and  our  knowledge  of  propa¬ 
gating  waves,  i.e.  this  representation  does  not  match  the  sign  convention  normally  used  to  represent 
the  direction  of  a  propagating  plane  wave. 

To  match  this  convention,  the  Fourier  transform  is  redefined.  The  4-D  case  (using  the  prop¬ 
agating  plane  wave  as  the  argument)  becomes 

/oo  /•oo 

/  5(x,t)e-J(‘"<-*'-*)dxdt.  (2.8) 

-00  — 00 

Notice  the  sign  is  only  changed  on  the  position  term,  x.  The  equation  now  correctly  reflects  the 
direction  of  the  plane  wave  from  basic  electromagnetics. 

2,2,1  Wavenumber-Frequency  Space,  The  (k,a>)  space  obtained  from  the  Fourier  trans¬ 
form  is  referred  to  as  wavenumber-frequency  space.  The  vector  k  is  the  wavenumber  vector.  This 
vector  represents  the  spatial  frequency  of  the  wave,  while  lj  represents  the  familiar  temporal  fre¬ 
quency.  Applying  the  4-D  Fourier  transform  to  the  propagating  plane  wave  of  Eqn,  (2.1)  uncovers 
an  important  property  of  the  frequency  spectrum  of  plane  waves.  The  sequence  of  operations 
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uncovering  this  property  is 


5(k,w) 


giK<-ko.x)g-j(a;t-k.x)^jj 


-Wo)e-^^^<’-^>^dx 


g-j(kc-k).x^3j 


^(a;  -  Wo)^(ko  -  k). 


(2.9) 


This  is  a  four  dimensional  impulse  at  k  =  ko  and  u  =  Uo»  Therefore,  every  point  in  (k,c«;)  space 
corresponds  to  a  plane  wave  in  (x,t)  space  with  a  particular  direction  and  frequency.  For  this 
reason,  antenna  beamforming  lends  itself  directly  to  multidimensional  signal  processing.  Designing 
a  bandpass  filter  in  wavenumber-frequency  space  transforms  back  to  space-time  as  an  antenna  beam 
passing  only  the  propagating  waves  with  a  particular  direction  and  frequency. 


2.2.2  Space-Time  Signal  Decomposition.  Another  important  concept  involves  the  de¬ 
composition  of  a  signal  in  (x,t)  space.  Any  space-time  signal  5(x,t)  can  be  decomposed  into  a 
superposition  of  propagating  plane  waves.  This  is  shown  by  the  4-D  inverse  Fourier  transform, 


=  7^  r  r  S{k,u;)^^^;;^dkdw, 


Plane  Wave 


(2.10) 


where  5(k,  cj)  is  the  4-D  Fourier  transform  of  the  space-  time  signal  s{x,  t).  Notice  the  exponential  in 
the  expression  exactly  resembles  the  expression  for  the  propagating  plane  wave  given  in  Eqn.  (2.1). 
^(kja;)  is  viewed  as  the  complex  magnitude  of  each  individual  plane  wave.  The  integrals  basically 
sum  all  of  the  individual  propagating  plane  waves  in  each  direction  to  form  the  space-time  signal. 
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Figure  2.1  Wavenumber-Frequency 
space  representation  of  space-time 
signals  with  the  same  frequency. 


Figure  2.2  Wavenumber-Frequency 
space  representation  of  space-time 
signals  with  the  same  velocity  of 
propagation. 


different  direction 


Figure  2.3  Wavenumber-Frequency 
space  representation  of  space-time 
signals  with  the  same  direction  of 
propagation. 


2.2.3  Wavenumber-Frequency  Space  Translation.  Using  the  fundamental  ideas  just  ex¬ 
plained,  the  plots  shown  in  Figs.  2.1  through  2.3  help  illustrate  how  the  (k,a;)  space  relates  back 
to  (x,  t)  space.  These  figures  ignore  spatial  frequency  in  the  z  direction  (kz  from  the  wavenumber 
vector)  in  order  to  allow  representation  in  three  dimensions. 

To  see  how  these  plots  are  related,  examine  the  relationship  between  velocity  of  propagation, 
temporal  frequency,  and  spatial  frequency.  Looking  at  the  units  is  the  easiest  way  to  understand 
the  relationship  between  these  three  quantities, 


(2.11) 
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Holding  the  temporal  frequency  constant  without  any  constraint  on  the  velocity  of  propaga¬ 
tion  or  direction  of  the  wave  results  in  a  horizontal  plane  in  wavenumber-frequency  space.  This 
is  shown  in  Fig.  2.1.  Forcing  the  velocity  of  propagation  to  be  constant  as  in  Fig.  2.2  results  in 
a  cone.  As  the  temporal  frequency  u  increases,  the  magnitude  of  the  spatial  frequency  must  also 
increase  at  the  same  rate  in  order  for  the  velocity  of  propagation  to  be  constant.  Hence,  the  area 
of  wavenumber- frequency  space  described  by  this  relationship  is  the  surface  of  a  cone. 

The  final  relationship  of  interest  is  signals  propagating  in  the  same  direction.  These  translate 
to  the  space  shown  in  Fig.  2.3.  The  direction  of  the  signal  is  described  by  the  wavenumber/spatial 
frequency  vector  k.  The  direction  of  the  vector  is  the  same  as  long  as  it’s  individual  components  are 
all  operated  on  in  a  like  manner.  Any  increase  in  kx  must  have  a  corresponding  increase  in  ky  (the 
ratio  remains  constant).  Therefore,  one  direction  of  propagation  is  described  by  a  vertical  plane 
(no  temporal  frequency  restriction)  in  wavenumber-frequency  space.  The  plot  shown  illustrates  two 
different  directions. 

The  intersections  of  these  plots  give  a  visual  interpretation  to  the  signal  processing  aspect  of 
beamforming.  The  object  of  beamforming  is  to  eliminate  undesired  signals  and  focus  only  on  signals 
of  interest.  Generally  the  desired  signal  is  one  of  a  particular  direction,  velocity  of  propagation, 
and/or  frequency.  Signals  with  the  same  velocity  of  propagation  and  direction  lie  along  a  line 
formed  by  the  cone  and  vertical  planes  shown  in  Fig.  2.2  and  2.3.  This  line  can  be  interpreted  as 
a  bandpass  filter  as  shown  in  Fig.  2.4.  By  also  choosing  a  frequency,  this  intersection  line  collapses 
down  to  a  point.  In  radar,  frequency  is  generally  of  no  concern  since  a  typical  radar  system  is 
trying  to  receive  the  ’’echos”  from  its  own  transmission.  The  existing  roll-off  of  the  amplifiers  is 
usually  adequate  to  remove  undesired  frequency  components. 

In  summary,  let  us  translate  some  of  the  basic  ideas  to  practical  concepts.  The  position 
vector  X  represents  a  three  dimensional  array  of  antennas  used  to  spatially  sample  the  signal.  The 
signal  applied  to  each  antenna  element  (used  synonymously  with  receiver)  is  the  spatially  sampled 
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Figure  2.4  The  intersection  of  the  constant  velocity  of  propagation  cone  and  the 
directional  plane  in  wavenumber-frequency  space. 

signal  5(xi,i),  where  i  denotes  the  antenna  element.  Temporal  (time)  sampling  provides  the  final 
dimension.  Superposition  allows  the  space-time  signal  to  be  broken  down  into  a  sum  of  individual 
plane  wave  components,  each  of  a  particular  frequency  and  direction.  The  appropriate  bandpass 
filter  removes  undesired  signal  components. 

2,S  Weighted  Delay  and  Sum  Beamformer 

The  objective  of  beamforming  is  to  isolate  signals  propagating  in  the  same  direction.  Fig¬ 
ure  2.5  graphically  illustrates  the  idea  of  forming  a  beam  or  passband  in  the  direction  of  the  target. 
Undesired  signals  are  attenuated  by  the  roll-off  of  the  antenna  pattern.  Ideally,  the  beamformer 
passes  the  desired  target  (direction  and  speed  represented  by  the  slowness  vector  ao)  and  totally 
rejects  the  undesired  signals. 


Undesired  Signal 


Figure  2.5  A  block  diagram  of  the  signal  environment. 


Assuming  all  waves  propagate  at  the  same  speed,  the  beamformer  becomes  the  intersection 
of  the  cone  in  Fig.  2.2  and  the  vertical  plane  in  Fig.  2.3.  This  describes  signals  of  any  frequency 
constrained  to  a  particular  direction  and  velocity  of  propagation.  The  assumption  of  a  constant 
velocity  of  propagation  Vp  is  valid  in  the  discussion  of  radar  since  the  wave  is  propagating  in 
free  space,  generally  assumed  a  homogeneous  medium.  The  wavenumber-frequency  space  diagram 
depicting  this  is  shown  in  Fig.  2.4.  The  dashed  lines  on  each  side  of  the  intersection  line  represent 
the  passband  of  the  beamforming  filter. 

As  mentioned  at  the  end  of  the  previous  section,  the  directional  vector  x  is  interpreted  as  the 
location  of  the  antenna  elements.  For  an  antenna  consisting  of  N  elements,  each  element  position 
is  described  by  a  separate  vector  Xi.  The  diagram  in  Fig,  2.6  shows  the  antenna  elements  lined  up 
in  the  y  direction,  the  convention  to  be  used  in  this  chapter.  Assuming  these  elements  are  spaced 
equal  distances  apart,  the  position  vector  describing  each  element  location  becomes 

Xi  =  idy  for  z  =  0  :  AT  -  1,  (2.12) 
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bf(f) 


Figure  2.6  Block  diagram  of  the  weighted  delay  and  sum  beamformer.  The  out¬ 
put  of  the  beamformer  as  a  function  of  time  is  designated  bf(t). 


where  d  is  the  interelement  spacing.  The  reference  element  is  assigned  the  zero  position  in  the 
coordinate  system,  referred  to  as  element  zero.  Using  this  notation,  the  received  signal  (denoted 
ri{t))  from  each  element  becomes  5(xi,t),  the  space-time  signal  sampled  at  position  Xi. 

The  block  diagram  in  Fig.  2.6  also  shows  a  weight  applied  to  each  received  signal.  At  this 
point  in  the  derivation  the  weight  is  taken  to  be  a  real  valued  constant  for  each  element,  although 
not  constrained  to  the  same  value  for  every  element.  Complex  weights  are  time  delays  and  hence 
are  constrained  to  the  delay  blocks.  Without  any  delays  on  the  elements,  the  beam  is  perpendicular 
to  the  face  of  the  antenna  array  (the  x  direction).  However,  placing  individually  selectable  complex 
exponential  weights  (time  delays)  after  each  antenna  element  steers  the  beam.  The  sum  and  1/N 
blocks  after  the  weight  and  delay  blocks  serve  to  average  the  signal  received  from  each  element. 
Called  the  weighted  delay  and  sum  beamformer,  this  beamformer’s  name  is  a  direct  derivative  of 
its  elements  of  operation. 


2.3J  Beam  Steering  Through  Time  Delays-  The  delay  concept  is  simple.  A  plane  wave 
propagating  in  some  direction  other  than  perpendicular  to  the  array  arrives  at  different  times 
on  each  element.  Applying  the  appropriate  delay  in  each  element’s  path  before  summing  all  the 
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signals  together  allows  a  reference  point  on  the  signal  to  be  in  phase  at  the  sum  block.  Hence, 
it  adds  constructively  and  produces  a  large  output  from  the  beamformer.  However,  a  plane  wave 
propagating  in  a  different  direction  experiences  phase  differences  when  it  reaches  the  sum  block.  A 
weak  beamformer  output  signal  bf(t)  results  from  the  destructive  addition.  Delays,  given  by 

ri  =  -ao-Xi  (2.13) 


pass  waves  in  the  direction  of  the  slowness  vector  [5,  page  294]. 

For  example,  examine  the  diagram  in  Fig.  2.7.  Drawing  a  reference  line  perpendicular  to 
the  direction  of  the  propagating  wave  and  based  at  element  zero  shows  the  increase  in  distance  a 
reference  point  on  the  wave  travels  to  reach  element  one.  The  time  to  travel  this  distance  is  the 
negative  delay  required  in  element  one.  Calculating  this  delay  is  simply  a  matter  of  dividing  the 
distance  traveled  by  the  velocity  of  propagation.  In  this  simple  example,  all  of  the  elements  are 
lined  up  vertically  (y  direction)  and  equally  spaced.  Therefore,  once  the  required  delay  for  element 
one  is  found  the  rest  of  the  delays  are  integer  multiples  of  it. 

Given  d  as  the  distance  between  the  elements,  the  additional  travel  distance  a  (see  Fig.  2.7) 
is 

a  =  dcos(0).  (2.14) 

The  velocity  of  propagation  is  the  temporal  frequency  of  the  wave  multiplied  by  the  wavelength. 


c  =  /oA. 


(2.15) 


Desired  Incoming  Wave 


Figure  2.7  Calculation  of  the  delay  required  in  each  element  to  steer  the  beam. 

Therefore,  the  delay  required  in  element  one  for  the  reference  point  to  arrive  simultaneously  at  the 
sum  block  from  elements  zero  and  one  is 


__  a 

c 

_  d  cos{9) 

_  27r  d  cos(0) 

A  27r/o 

_  dcos{9)  27r 


Spatial  Frequency 
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Recognizing  the  spatial  frequency  term,  the  equation  becomes 


T  =  _ll^dcos(0) 

LOq 

=  -|Q!o|dcOS(0) 

=  -ao-Xi,  (2.16) 

where  the  final  step  shows  the  dot  product  between  the  two  vectors  (the  same  as  a  vector  multi¬ 
plication)  and  gives  the  expression  provided  by  Dudgeon  and  Mersereau  [5,  page  294],  Examining 
the  element  position  vector  (z  =  1  for  element  one)  in  this  example, 

Xi  =  [0  id  of,  (2.17) 

shows  how  the  dot  product  of  Eqn.  (2.16)  collapses  to  the  expression  on  the  previous  line. 

Obviously,  a  negative  time  delay  is  actually  an  advance  and  hence  not  possible  in  the  real 
world.  Adding  the  absolute  value  of  rjv-i  to  all  of  the  delay  blocks  produces  the  more  conventional 
notion  of  a  positive  time  delay.  In  this  case,  r^-i  becomes  the  delay  block  equal  to  zero  instead  of 
To. 

2.3.2  Antenna  Array  Pattern.  Examination  of  the  weighted  delay  and  sum  beamformer 
block  diagram  in  Fig.  2.6  allows  a  derivation  of  a  mathematical  expression  for  the  output.  The 
beamformer  output  is  expressed  as 

.  N-l 

2=0 

where  ri{t  —  ri)  is  the  received  signal  sampled  at  a  particular  position  x^.  Assuming  the  incoming 
signal  is  a  single  plane  wave,  then  (using  the  plane  wave  expression  from  Eqn.  (2.5)  and  the  delay 
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expression  from  Eqn.  (2.13))  the  beamformer  output  becomes 


m 


1 

—  ^  Wia(xj,i  +  ao-Xi) 
i=0 

1 

—  yj.QMt+ac-Xi-a-Xi) 

^  i=0 
1  ■ 

i=0 

1 

i=0 


W(k-w«„) 


(2.19) 


where,  the  underbraced  term  is  defined  as  W(k  —  loao),  the  antenna  array  patem  [5]. 

The  above  development  of  the  antenna  array  pattern  assumed  the  space-time  signal  consisted 
of  only  one  propagating  wave.  However,  we  have  already  shown  the  space-time  signal  is  composed 
of  many  plane  waves  (see  Eqn.  (2.10)).  Using  this  more  general  case,  the  received  signal  becomes 
Eqn.  (2.10)  evaluated  at  Xj  and  the  necessary  beam  steering  delay  from  Eqn.  (2.13), 


ri{t  -  Tj) 


s{Xi,t-Ti) 


1 

(27r)4 


5(k,  dw. 


(2.20) 


Substituting  this  equation  into  the  generic  expression  for  the  delay  and  sum  beamformer  (Eqn.  (2.18)) 
results  in 


1  (  1  f°°  r°°  1 


(2.21) 
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Rearranging  the  order  of  the  summation  and  integral  highlights  the  expression  for  the  antenna 
array  pattern, 


1  ^00  poo  (  1  ^”1  'I 

-j  poo  pOO 

=  -zr-u  /  /  5(k,w)W(k-a;ao)e^"*dkda;. 

^I'K)  J  J  —  OO 


(2.22) 


Examine  the  response  of  the  weighted  delay  and  sum  beamformer  for  a  space-time  signal  with 
all  components  propagating  in  the  same  direction.  The  space-time  signal  becomes 


s{x,t)  =  v{t  —  a  •  x) 


(2.23) 


for  the  signal  v{t)j  where  the  component  -a  •  x  is  the  time  delay  from  Eqn.  (2.13).  As  before,  a 
represents  the  direction  of  the  incoming  signal  and  ao  represents  the  steered  direction  of  the  beam. 
The  Fourier  transform  of  the  incoming  signal  (as  defined  to  correctly  reflect  the  direction  of  a  plane 
wave)  involves  using  the  time  shift  property. 


5(k,a;)  =  /  /  3{x,t)e  dxdt 

J  —  00  «/  —  00 

/oo  poo 

/  v(t  —  a  •  dxdt 

-00  J  —  oo 

/oo  poo 

/  v{t  —  a-  dtdx 

-00  J  —  00 

=  V{u)  f 

J  —  00 

=  y(«)%a-k). 


(2.24) 


This  shows  the  (k,  u)  space  of  a  propagating  wave  is  a  single  point. 

Expressing  the  beamformer  output  as  a  function  of  the  space-time  signal  with  all  components 
propagating  in  the  same  direction  shows  the  effect  of  the  antenna  array  pattern.  The  sifting  property 
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of  the  delta  function  greatly  simplifies  the  mathematics  involved  in  reducing  the  expression, 


bf(f)  =  r  r  V{u)6{wa-k)Wik-(jjao)e^‘^*dkcLj 

(zTT)^  J-oo  J-oo 
1  poo 

=  /  V{u)W{ij{a-ao))e^‘^^(Lj. 

"TT  J—oo 


If  the  steered  beam  direction  is  the  same  as  the  direction  of  the  incoming  signal,  a  =  a©,  then  the 
above  integral  collapses  to  the  inverse  Fourier  transform  of  the  original  signal,  F(a;),  multiplied  by 
a  constant, 


J — oo 

=  Ty(0)u(t). 


(2.26) 


As  expected,  the  signal  is  undistorted.  If  a  7^  aoj  then  the  antenna  array  pattern  attenuates  the 
signal.  Therefore,  --  uao)  can  be  viewed  as  the  antenna  pattern. 

2.3.3  Linear  Antenna  Array  Example.  As  an  example  of  the  computation  of  an  antenna 
pattern  for  a  weighted  delay  and  sum  beamformer,  consider  a  linear  antenna  array  with  N  elements 
aligned  vertically  (y  direction).  Given  the  weights  w  are  uniform  and  real  valued  for  all  of  the 
antenna  elements,  calculate  the  antenna  pattern.  Assume  the  antenna  elements  are  equally  spaced 
a  distance  d  apart  and  the  beam  is  steered  toward  the  direction  k©. 
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First,  we’ll  determine  the  delays  necessary  at  each  element  to  steer  the  beam  in  this  direction. 


Referencing  Eqn.  (2.13)  gives 


T'i  —  Q^q  *  Xf 
=  “O^oko  •  Xi 

—  [  ^oa;  ^oy  ^02]  [0  0] 

=  —idwokoy  =  —idaoy  for  i  =  0  :  iV  —  1. 


Finding  a  closed  expression  for  the  antenna  array  pattern  is  merely  a  matter  of  manipulation. 
The  antenna  pattern  for  the  steered  antenna  beam  given  a;  =  a;^  such  that  ujoto  =  ij^oOto  =  k©  is 

W{\i 


Substituting  the  values  in  the  exponential  argument  and  simplifying  results  in 


W{k-ko)  = 


^  r,-j{ky-kov)id 

N  ^ 

t=l 

^  1  —  (e-j(*Jy”*ov)^)^ 

jV'  ][  —  Q—j{ky—koy)d 

'll)  1  ^  ^—j(ky -‘koy)dN 

N  1  —  Q~^j{ky—koy)d 

*11)  g  J(,ky  koy')dN ^ji^ky  koy)dN/2 _ g  ji.ky  koy')dN/2 

jV  Q—j{ky—koy)d/2  Qj{ky—koy)d/2  —  Q—j{ky—koy)dl2 

w  sin((fcy  —  koy)dN/2) 

N  sin((A:j,  -  koy)d/2) 


A  plot  of  the  magnitude  for  this  beamformer  (with  w  —  1  and  ko  =  0)  is  shown  in  Fig.  2.8. 
Examining  the  arguments  of  the  sine  terms  shows  the  pattern  is  periodic  in  2n/d  radians/meter. 
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Antenna  Array  Pattern 


Figure  2.8  Antenna  pattern  magnitude  for  the  weighted  delay  and  sum  beam- 
former,  10  elements  spaced  0.1m  apart,  and  no  beam  steering. 

The  first  null  occurs  at  27r/(Ard).  The  repeated  main  lobes  are  called  grating  lobes.  Figure  2.9 
illustrates  steering  the  beam  to  the  koy  ==  it /d  direction. 

2.3.4  Beam  Shaping  Through  Windowing.  The  weights  in  the  beamformer  perform  the 
function  of  a  multidimensional  filter  window.  In  the  example  above,  the  constant  valued  weights 
correspond  to  a  one  dimensional  rectangular  window,  explained  pictorially  by  Fig.  2.10.  The  window 
affects  the  width  of  the  mainlobe  and  the  level  of  the  sidelobes.  Ideally,  the  antenna  pattern  should 
be  a  delta  function,  i.e.  a  mainlobe  of  infinitesimal  width  and  no  sidelobes.  This  translates  to  the 
ability  to  pick  out  targets  no  matter  how  close  together  they  are.  However,  this  is  not  physically 
realizable.  Using  the  rectangular  window  in  the  previous  example  produced  a  reasonably  narrow 
mainlobe  resulting  in  the  antenna  receiving  a  small  range  of  directions  at  the  expense  of  higher 
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One  Dimensional  Window  Coefficients 


Figure  2.10  Pictorial  representation  of  one  dimensional  rectangular  and  Ham- 
ming  windows  used  for  the  receiver  weights. 

sidelobes  (reference  Fig.  2.8).  The  advantage  of  the  rectangular  window  is  simplicity  through 
weights  of  constant  value,  shown  by  the  relatively  simple  expression  developed  for  the  antenna 
pattern  in  the  previous  example. 

Other  windowing  options  include  the  Hanning,  Hamming,  Bartlett,  Kaiser,  Chebyshev,  and 
Gaussian  windows.  The  expression  for  the  weights  using  the  one  dimensional  Hamming  window  [16, 
page  183]  is 


Wi  =  0.54  +  (1  -  0.54)  cos 


(2.27) 
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Antenna  Array  Pattern 


Figure  2.11  Hamming  window  applied  to  the  weights  of  a  linear  antenna  array 
of  10  elements  space  0.1m  apart. 


The  shift  on  the  cosine  term  is  added  to  center  the  window  on  the  receiver  array.  The  derivation 
for  the  previous  linear  antenna  array  example  using  a  Hamming  window  for  the  weights  is 


W{k-ko) 


i=0 


Due  to  the  more  complicated  element  weighting,  the  equation  cannot  be  simplified  as  in  the 
rectangular  window  case.  Figure  2.11  shows  the  plot  of  this  expression.  Comparing  this  with 
Fig.  2.8  highlights  the  trade-off  between  a  narrow  mainlobe  and  low  sidelobes.  Also  notice  the  loss 
in  mainbeam  gain. 
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Two  dimensional  windows  can  be  applied  to  the  weights  to  control  the  antenna  pattern  of  a 
rectangular  antenna  array.  There  are  also  windowing  techniques  applicable  to  a  hexagonally  spaced 
array,  corresponding  to  a  hexagonal  sampling  lattice  in  the  spatial  domain.  As  long  as  the  antenna 
array  is  on  a  periodic  grid,  the  design  of  the  antenna  pattern  using  multidimensional  filters  is  rather 
straightforward.  However,  nonuniform  sampling  significantly  increases  the  difficulty  in  designing 
the  receiver  weighting  function. 

Summary 

In  this  chapter,  the  fundamental  concepts  of  beamforming  were  introduced.  The  introduction 
of  the  plane  wave  model  and  antenna  beamformer  development  showed  the  relevance  of  signal 
processing  techniques.  The  wavenumber-frequency  plots  illustrated  the  fundamental  relationships 
between  the  wavenumber-frequency  and  space-time  domains. 

From  these  basic  ideas,  the  weighted  delay  and  sum  beamformer  was  developed  and  exam¬ 
ined.  This  beamformer  serves  as  an  excellent  introduction  to  antenna  arrays.  Using  a  linear  array 
example,  the  antenna  pattern  was  derived  and  examined.  The  receiver  weights  were  designed  using 
windowing/signal  processing  techniques.  The  use  of  a  different  window  was  shown  to  reduce  the 
sidelobes  of  the  antenna  pattern,  at  the  expense  of  mainbeam  width  and  gain.  It  should  be  noted 
the  weighted  delay  and  sum  beamformer  can  be  further  generalized  by  allowing  the  filter  weights 
to  depend  on  frequency.  This  configuration  is  referred  to  as  the  filter  and  sum  beamformer. 

All  of  the  equations  to  this  point  treat  time  as  a  continuous  variable.  Implementation  of  a 
beamforming  system  with  digital  hardware  necessitates  sampling  the  signal  temporally.  One  issue 
brought  forth  by  the  sampling  operation  is  a  distortion  of  the  array  pattern  when  the  time  delays 
used  to  steer  the  beam  are  not  multiples  of  the  sampling  period.  Dudgeon  and  Merserau  [5,  pages 
303-307]  discuss  this  in  more  detail. 
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III.  Space-Time  Adaptive  Processing  Primer 

3.1  Introduction 

Airborne  radars  hold  many  advantages  over  ground  based  radars,  however  they  also  suffer 
several  disadvantages.  Higher  clutter  levels  due  to  steep  aspect  angles  to  the  ground  make  it  difficult 
to  detect  and  track  low  flying  targets.  Airborne  radars  are  also  hindered  by  the  constant  motion  of 
the  airborne  platform  relative  to  the  ground  and  the  changing  terrain,  adding  a  variable  nature  to 
the  clutter  returns.  The  typical  airborne  radar  environment  consists  of  spectrally  variable  ground 
clutter  in  both  the  mainlobe  and  sidelobes  of  the  radar  and  also  interference  due  to  jammers.  The 
aircraft  velocity  adds  Doppler  spreading  to  returns  processed  by  the  radar. 

Radar  engineers  have  made  numerous  attempts  to  improve  the  detection  capabilities  of  air¬ 
borne  radars.  One  of  the  most  common  is  to  reduce  the  effects  of  clutter  through  the  use  of  lower 
sidelobes.  However,  practical  limitations  of  most  antenna  designs  limit  this  to  the  order  of  only 
-30  to  -40  dB,  allowing  Doppler  spread  sidelobe  clutter  returns  to  be  strong  enough  to  completely 
mask  low-observable  targets  [12].  Although  antennas  with  lower  sidelobe  levels  are  achievable,  the 
costs  associated  become  prohibitive.  Also,  lower  sidelobes  do  not  offer  any  gains  in  the  ability  to 
detect  a  target  amid  jamming. 

Another  option  is  non-adaptive  filtering.  In  a  ground  based  radar,  reasonably  constant  clutter 
patterns  allow  non-adaptive  filtering.  However,  the  Doppler  spread  and  variable  nature  of  the 
clutter  in  airborne  radars  make  this  option  better  suited  for  other  applications. 

In  the  last  decade,  improvements  in  digital  technology  make  adaptive  filtering  algorithms  very 
attractive.  An  adaptive  filtering  algorithm  can  have  significant  impact  on  radar  performance  both 
amid  clutter  and  jamming.  Not  only  is  target  acquisition  improved,  but  the  overall  design  of  the 
radar  is  affected.  For  example,  the  ability  to  filter  sidelobe  clutter  reduces  the  need  for  extremely 
low  antenna  sidelobes.  Another  example  is  the  Pulse  Repetition  Frequency  (PRF)  design  trade-off. 
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Typically  the  PRF  is  chosen  to  help  eliminate  the  effects  of  clutter,  if  clutter  can  be  minimized 
through  the  use  of  filtering,  the  PRF  can  be  chosen  to  give  the  radar  more  unambiguous  range. 

Space-Time  Adaptive  Processing  (STAP)  refers  to  this  class  of  adaptive  filtering  algorithms. 
Using  digital  beamforming  concepts,  STAP  is  defined  by  Ward  [21]  as  the  simultaneous  combination 
of  signals  from  the  elements  of  an  array  antenna  (the  spatial  domain)  and  the  multiple  pulses  of 
a  coherent  radar  waveform  (the  temporal  domain)  to  suppress  interference  and  provide  target 
detection.  STAP  is  inherently  an  extension  of  digital  beamforming  architecture  and  principles. 

Chapter  II  introduced  the  basics  of  beamforming  and  serves  as  the  foundation  for  the  prin¬ 
ciples  encountered  in  this  next  development.  As  explained,  the  beamforming  operation  occurred 
through  the  use  of  weights  on  each  antenna  element.  The  output  of  the  summed  and  weighted  ele¬ 
ments  resulted  in  a  single  data  stream  containing  signals  from  propagating  waves  in  the  direction  of 
interest.  STAP  differs  from  basic  antenna  beamforming  in  that  no  weights  are  used  on  the  elements 
themselves.  The  received  space- time  signal  from  each  element  is  stored  for  later  processing.  The 
weighting  functions  are  applied  as  part  of  the  processing  algorithm. 

The  development  begins  by  defining  the  radar  geometry.  Models  are  introduced  for  noise, 
jamming,  and  clutter  along  with  the  target  return  signal.  The  nomenclature  and  architecture  of 
STAP  is  introduced  as  the  development  progresses.  Chapter  II  concepts  are  tied  together  to  form 
these  models.  Ward’s  [21]  work  is  the  basis  for  the  developments  in  this  chapter. 

3,2  The  Physical  and  Electromagnetic  Environment 

3.2,1  Physical  Geometry.  The  radar  antenna  array  for  all  examples  in  this  chapter  is 
assumed  linear  with  uniformly  spaced  elements  a  distance  d  apart,  each  with  their  own  receiver. 
The  orientation  of  the  array  is  along  the  x-axis,  as  shown  in  Fig.  3.1a.  Therefore,  the  vector  to  the 
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Figure  3.1  (a)  Uniformly  spaced  antenna  element  array  geometry  and  (b)  radar 

platform  geometry. 


first  antenna  element  is 


d  =  dx, 


(3.1) 


where  x  is  the  unit  vector  in  the  x  axis  direction.  Using  this  expression,  the  vector  to  the 
element  is 


Tn  =  nd. 


(3.2) 


As  used  in  Fig.  3.1b,  the  vector  Vo  shows  the  magnitude  and  direction  of  the  platform  velocity. 
This  figure  also  shows  the  corresponding  radar  platform  geometry. 

The  unit  vector  k(^,0)  shown  in  Fig.  3.1b  represents  the  direction  to  a  point  of  interest 
(or  target).  This  direction  corresponds  to  the  direction  of  k,  the  spatial  frequency  vector  from 
Chapter  II,  see  Eqn.  (2.1).  The  angular  variables  <f>  and  6  respectively  represent  the  azimuth  and 
elevation  to  the  point  of  interest  (target,  clutter,  jammer,  etc.).  The  relationship  of  Ward^s  [21] 
coordinate  system  to  the  Cartesian  coordinate  system  is  described  by 

k(0, 9)  =  cos  9  sin  0x  +  cos  9  cos  (f>y  +  sin  0z,  (3.3) 
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where  y  and  z  again  represent  unit  vectors  along  the  Cartesian  coordinate  system  y  and  z  axes. 


3.2,2  Signal  and  Radar  Description.  Now  that  the  physical  geometry  of  the  radar  is 
defined,  let  us  examine  the  electromagnetic  properties  of  the  radar.  The  radar  is  assumed  to  be 
Pulse  Doppler  (PD),  therefore  target  velocity  information  is  available  to  the  processor. 

The  sinusoidal  waveform  transmitted  by  the  radar  is  modeled  mathematically  by 

=  (3.4) 

where  the  complex  exponential  and  random  phase  shift  represent  the  sinusoidal  carrier,  at  the 
transmitted  pulse  amplitude,  and  u{t)  the  envelope  function.  The  envelope  function  defines  the 
pulse  width  and  Pulse  Repetition  Interval  (PRI)  characteristics.  A  series  of  individual  pulses 
summed  to  form  the  pulse  train  of  the  radar  defines  the  envelope  function, 

M-l 

u{t)  =  ^2  Up{t  —  mTr).  (3.5) 

m=0 

The  fundamental  pulse  is  time  shifted  by  multiples  of  the  PRI,  Tr.  The  summation  of  a  finite 
number  of  terms,  M,  shows  the  Coherent  Processing  Interval,  or  CPI.  This  expression  implicitly 
defines  M  as  the  number  of  pulses  in  the  CPI. 

One  should  be  careful  to  distinguish  between  the  coherent  processing  interval  and  the  coherent 
integration  time.  We  will  find  later  that  the  output  of  the  STAP  algorithm  is  described  by  2/  = 
where  w  is  the  adaptive  weight  vector  and  x  is  the  incoming  space-time  snapshot  (discussed  later 
in  detail).  This  expression  indicates  the  coherent  integration  time  since  the  inner  product  of  the 
weight  vector  and  space-time  snapshot  is  simply  the  weighted  sum  of  the  incoming  data  samples, 
a  form  of  integration. 
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The  CPI  defines  how  much  data  the  adaptive  processor  uses  in  determining  the  adapted 
pattern.  Obviously,  a  smaller  CPI  allows  rapid  modifications  to  the  antenna  pattern  in  response 
to  the  changing  environment.  However,  making  it  too  small  does  not  give  the  processor  enough 
information  about  the  statistics  of  the  environment,  producing  poorly  formed  antenna  patterns. 
Making  the  CPI  too  large  results  in  slow  response.  This  is  one  of  the  primary  design  considerations 
in  a  STAP  system. 

The  transmit  signal  amplitude  at  is  defined  to  capture  the  energy  in  the  pulse.  Therefore, 
the  fundamental  pulse  energy  Eu  is  normalized  to  one, 

Eu=  =  (3.6) 

Jo 

Consequently,  the  total  energy  transmitted  in  a  single  pulse  Ep  meets  the  definition  for  at, 

r  mfdt  =:al  (3.7) 

Jo 

From  this  expression,  the  energy  within  one  CPI  is  easily  determined  by  integrating  over  the 
CPI  length.  Since  the  focus  of  STAP  algorithms  is  the  response  and  pattern  after  one  CPI,  it  is 
customary  to  consider  only  the  signal  transmitted  within  one  CPI  as  being  the  entire  transmitted 
signal.  Hence,  the  energy  within  one  CPI  is  synonymous  with  the  transmitted  signal  energy  Ef, 

tMTr 

Et=  \s{t)fdt  =  Mai  =  MEp.  (3.8) 

Jo 

The  signal  received  by  each  individual  antenna  element  is  the  transmitted  signal  (Eqn.  (3.4)) 
with  a  time  delay  from  the  target  to  the  element,  Tn,  and  a  Doppler  shift  ft  due  to  the  target’s 
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relative  velocity  between  the  target  and  the  platform,  ft  =  2'yt/Ao,^ 

S„{t)  =  s{t  -  (3.9) 

Inserting  the  transmitted  pulse  expression  gives 

Sn{t)  =  aMt  -  (3.IO) 

where  ar  represents  the  received  pulse  amplitude.  This  parameter  is  introduced  because  the  atmo¬ 
sphere  and  Radar  Cross  Section  (RCS)  of  the  target  attenuate  the  amplitude  of  the  pulse,  hence  it 
is  typically  much  less  than  Of. 

If  the  Doppler  shift  is  greater  than  the  Pulse  Repetition  Frequency  (PRF),  ft  >  fr,  then  the 
velocity  of  the  target  is  ambiguous.  This  is  shown  in  Fig.  3.2.  Once  the  target  Doppler  shift  has 
increased  past  the  first  PRF  line  located  at  /r,  the  radar  can  no  longer  unambiguously  determine 
the  target  velocity.  For  this  reason,  target  Doppler  frequency  is  normalized  to  the  PRF, 

=  ^  =  ftTr.  (3.11) 

Jr 

Using  the  normalized  Doppler  frequency  u),  a  quick  glance  determines  if  the  Doppler  shift  is  un¬ 
ambiguous  (normalized  Doppler  less  than  one)  or  ambiguous  (normalized  Doppler  greater  than 
one). 

Similar  to  Chapter  II,  the  signal  delay  is  expressed  as  a  delay  from  the  antenna  element  of 
interest  to  the  reference  element  (the  first  element  to  receive  the  return  signal),  r^,  and  a  delay 

^The  Doppler  shift  is  assumed  equal  at  all  antenna  elements. 
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Time  Domain 


Pulse  Spectrum  and  PRF  Lines 


equal  to  the  round  trip  time  to  the  target,  rt, 


Tn  =  Tt  +  2t^. 


(3.12) 


Keeping  the  physical  geometry  of  the  scenario  under  consideration,  we  see  the  antenna  element’s 
delay  is  a  function  of  the  target’s  elevation  0*,  azimuth  and  range  Rt.  The  round  trip  time  to 
the  target  is  a  simple  function  of  the  velocity  of  propagation  c  and  range  Rt, 


n 


(3.13) 


The  differential  delay  from  the  reference  element  to  the  element  is 


^  ^  _nHeu'k)-i  ^ 


(3.14) 
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Figure  3.3  Antenna  array  element  block  diagram.  ADC  is  an  analog  to  digital 
converter. 


Related  to  this  expression  is  the  target  spatial  frequency, 


^  ^  cos  0t  sin  (^t. 


(3.15) 


The  above  equation  for  the  spatial  frequency  relates  directly  to  the  |ko|  found  in  Chapter  II 
Eqn.  (2.3)  for  an  array  oriented  along  the  x  axis. 

Using  these  expressions,  the  received  signal  shown  in  Eqn.  (3.10)  can  be  rearranged  to  show 
the  radar  ambiguity  function.  Since  each  antenna  element  has  its  own  receiver  and  matched  filter 
as  shown  in  Fig,  3.3,  the  matched  filter  output^  for  the  n**  element  signal  is  given  as  [21], 

X„(f)  =  ^  -  r*  -  mT„  Qt),  (3.16) 

m=0 

where  xi^t,  ft)  is  the  radar  ambiguity  function  defined  by  Skolnik  [17,  page  412]  as 

Ixin,  /t)P  =  I  r  ui0)u*{0  +  ' .  (3.17) 

I*/ — oo 


The  ambiguity  function  is  used  to  assess  the  properties  of  the  transmitted  waveform  in  terms 
of  its  target  resolution,  measurement  accuracy,  ambiguity,  and  response  to  clutter  [17].  Section  11,4 

^Ward’s  expression  assumes  rt  »  2r^. 
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of  Skolnik’s  text  [17]  is  an  excellent  introduction  to  the  usefulness  of  the  radar  ambiguity  function. 
For  the  purposes  of  this  discussion,  the  following  property  is  of  interest. 

The  ambiguity  function  evaluated  at  (0,0)  represents  twice  the  energy  in  the  return  signal. 
Because  of  the  waveform  normalization  shown  in  Eqn.  (3.6)  this  is  simply  one, 

|x(0,0)|2  =  (2E)2  =  1.  (3.18) 

Let  us  consider  the  return  from  a  target  at  a  particular  range  gate,  L  First,  discretize  time 
through  a  sampling  operation  at  each  PRL  In  other  words,  the  continuous  variable  t  becomes 

tm  mTr  for  m  =  0  :  M  —  1,  (3.19) 

where  the  round  trip  delay  r<  is  introduced  to  simplify  the  next  step.  The  sampled  target  return 
from  a  particular  range  gate  is 

M-l 

Xnitm)  =  Xnm  =  ^  ((m  “  q)Tr,  ft) ,  (3.20) 

q=0 

where  Qt  is  the  normalized  Doppler  frequency  and  the  range  gate  dependence,  Z,  has  been  suppressed 
from  the  subscript.  Prom  this  point  forward,  we  always  assume  the  item  of  interest  is  for  a  particular 
range  gate/cell  allowing  suppression  of  the  I  subscript.  It  is  assumed  no  target  range  ambiguities 
exist,  therefore  the  x  function  is  zero  everywhere  except  when  g  =  m.  This  allows  the  simplification 

ft).  (3.21) 
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Furthermore,  the  product  of  the  target  Doppler  shift  and  pulse  width  is  assumed  much  less  than 
one  allowing  the  approximation 

ftTp  <  1  implies  x(0,  /)  «  1.  (3.22) 

Using  this  approximation,  the  target  samples  at  the  antenna  element  and  pulse,  Xnm  can 
be  expressed  as 

^  (3  23) 

where  the  complex  random  amplitude  term  at  contains  the  phase  and  amplitude  terms  shown 
previously,  i.e.  at  =  are^'^ . 

5.5  Fundamentals  of  STAP 

As  mentioned  previously  in  Chapter  II,  Matlab®  notation  is  used  wherever  possible  to 
express  the  formulation  of  matrices  and  vectors.  In  particular,  X  =  [xi  X2]  places  the  two  column 
vectors  side  by  side  to  form  a  matrix,  and  x=vec(X)  places  the  columns  of  X  end-to-end  forming 
a  vector.  Lower  case  bold  letters  represent  vectors  and  upper  case  bold  represent  matrices. 

To  begin  the  fundamentals  of  STAP,  examine  Eqn.  (3.23).  First,  create  a  column  vector 
containing  each  antenna  element  sample  at  a  particular  range  gate  I  and  pulse  m  denoting  it  x^, 

[1;  ; . . .  ; 

=  (3.24) 

where  is  the  spatial  steering  vector  of  length  N  and  the  range  gate/cell  dependence  has  again 
been  suppressed  in  the  subscript.  Using  this  notation,  the  incoming  samples  can  be  written  as  the 
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matrix  X, 


X  = 


Xo  Xi  ...  XM-1 


(3.25) 


where  the  matrix  is  for  the  range  gate/cell  of  interest. 

Since  the  spatial  frequency  is  a  function  of  the  target  elevation  and  azimuth  angles  (see 
Eqn.  (3.15)),  the  spatial  steering  vector  can  be  written  as  an  explicit  function  of  these  variables, 


(3.26) 


Similarly,  the  temporal  steering  vector  bt{u)t)  is  defined  as  the  vector  of  exponentials  spanning 
the  M  pulses  in  the  CPI, 


b<(a>()  =  [1  ...  _ 


(3.27) 


The  relationship  between  these  two  vectors  forming  the  matrix  of  received  target  samples 
is  characterized  by  the  Kronecker  product,  denoted  ®.  See  Appendix  A  for  a  discussion  of  its 
properties.  The  Kronecker  product  of  the  temporal  and  spatial  steering  vectors. 


vt(cDi,i?<)  =  bi(a;f)  O  a-ti'&t),  (3.28) 

is  defined  as  the  space-time  steering  vector -vt,  where  the  above  equation  emphasizes  it’s  normalized 
Doppler  and  spatial  frequency  dependence.  Scaling  this  length  NM  vector  by  a  random  received 
amplitude  parameter  (subscripted  with  a  t  for  the  target)  gives 


Xt  = 


(3.29) 


3-11 


where  Xt  is  the  space-time  snapshot  for  the  target.  This  length  NM  column  vector  represents  the 
target  returns  received  by  each  antenna  element  {n  =  0  :  N  —  1)  due  to  each  transmitted  pulse 
(m  =  0  :  M  -  1)  at  a  particular  range  gate/cell  of  interest  L  Relating  this  Kronecker  product  to 
the  matrix  of  returns  shown  in  Eqn.  (3.25)  illustrates  it  is  simply  the  stacking  of  the  columns  of  X 
end  to  end  to  form  the  column  vector, 

X  =  vec(X).  (3.30) 

Since  the  target  is  not  the  only  echo  returned  to  the  radar,  an  undesired  component,  Xu  niust 
be  introduced  in  the  model. 


X  =  X<  +  X«-  (3.31) 

Our  model  will  consider  the  undesired  components  due  to  clutter,  noise,  and  barrage  noise  jamming. 

Let  us  give  a  quick  summary  of  the  important  points  of  the  development  so  far.  The  complex 
scalar  Xnm  contains  the  received  signal  from  a  particular  antenna  element  and  pulse  of  interest. 
Expanding  this  into  a  vector  across  all  of  the  antenna  elements  results  in  the  spatial  steering  vector 
a(i?)  containing  the  received  signals  from  every  antenna  element  at  the  pulse  of  interest.  Expanding 
Xnm  into  a  vector  across  all  pulses  results  in  the  temporal  steering  vector  b(a>)  containing  the 
received  signals  from  every  pulse  but  only  one  antenna  element.  The  space- time  snapshot  combines 
these  two  vectors  into  one  column  vector,  denoted  by  the  use  of  the  Kronecker  Product  notation. 
It  contains  all  of  the  received  signals  from  every  pulse  and  every  antenna  element  at  one  range  gate 
/.  Since  the  time  between  each  pulse  is  divided  into  L  range  bins,  there  are  L  total  space-time 
snapshots  in  the  CPI. 

Figure  3.4  is  a  visual  illustration  of  the  CPI  datacube.  The  shaded  cross  section  is  represen¬ 
tative  of  the  space-time  snapshot  where  the  columns  of  the  matrix/cross  section  form  the  snapshot 
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®  PiU  M-1 

Temporal  Samples 

Figure  3.4  The  CPI  datacube  of  dimension  M  x  N  x  L, 


vector.  As  mentioned  previously,  the  number  of  pulses  M  in  the  CPI  define  its  length.  Therefore, 
the  time  duration  of  the  CPI  is  MTr. 

Also,  notice  the  raw  data  is  stored  from  each  antenna  element.  No  weighting  is  applied  to  the 
received  signals  as  explained  in  Chapter  II.  The  weighting  is  performed  after  the  data  is  stored. 
This  is  the  function  of  the  space-time  adaptive  processor.  Withholding  the  weights  until  after  the 
data  in  the  CPI  is  known  allows  the  statistics  of  the  data  to  drive  their  calculation. 

The  next  step  in  the  development  of  STAP  fundamentals  is  analysis  of  the  undesired  com¬ 
ponent  Xu  Eqn.  (3.31).  This  component  consists  primarily  of  three  sources:  thermal  noise, 
jammers,  and  clutter.  The  following  sections  examine  each  of  these  in  detail. 

3.3.1  Thermal  Noise  Model.  This  model  begins  by  assuming  the  only  noise  present  is  due 
to  the  internally  generated  thermal  noise  of  the  receiver.  This  noise  is  taken  to  be  white,  therefore 
there  is  no  correlation  of  the  noise  between  pulses  and  it  is  also  uncorrelated  spatially  (between 
antenna  elements).  Let  Xnm  describe  the  noise  received  in  the  antenna  element.  Since  the  noise 
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is  assumed  white,  the  expected  value  of  received  noise  is 


£{Xn 


nm^n'm 


(3.32) 


where  represents  the  noise  power  in  a  single  element.  This  noise  power  is  given  by  the  noise  power 
spectral  density  No  multiplied  by  the  bandwidth  B  of  the  receiver  (since  this  analysis  assumes  a 
separate  receiver  on  each  antenna  element),  =  NqB,  The  Kronecker  delta  <5  in  the  above  equation 
is  defined  as 


1  n  =  n^ 

0  n  ^  7i'. 


(3.33) 


A  similar  expression  shows  the  correlation  of  the  noise  terms  from  pulse  to  pulse. 


£{x< 


nm*^nm* 


}  =  Cr^Sm-r, 


(3.34) 


If  we  let  Xn  represent  the  MN  xl  spatial  snapshot  of  noise  signals  on  each  antenna  element, 
the  noise  space-time  covariance  matrix  is  shown  by 

Rn  =  £{XnXn  }  =  ^  liVj  (3.35) 

where  I  is  the  identity  matrix  with  dimension  shown  by  the  subscript  and  is  the  previously  defined 
noise  power  per  antenna  element.  The  covariance  matrix  shown  above  serves  as  the  thermal  noise 
component  of  the  overall  covariance  matrix  R.  The  next  undesired  component  considered  is  the 
jamming  model. 

3,3.2  Jamming  Model  At  this  point,  only  noise  jamming  originating  from  a  long  range 
airborne  or  ground  based  platform  is  considered.  The  bandwidth  of  the  jammer  Bj  is  assumed 
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larger  than  the  radars  instantaneous  bandwidth  (the  inverse  of  the  time  domain  pulse  width  Tp), 


i.e. 


=  (3.36) 

In  order  to  counter  a  jammer  through  STAP,  its  signal  must  remain  correlated  across  the 
antenna  elements.  This  means  the  time  for  the  jamming  signal  to  travel  across  the  array  is  less 
than  the  inverse  of  its  bandwidth,  a  measure  of  the  correlation  time  for  the  jammer.  Figure  3.5 
helps  explain  this  point. 

Examining  the  inverse  Fourier  transform  of  the  jammer’s  power  spectral  density  envelope 
shows  an  increase  in  the  time  domain  width  as  the  jammer  bandwidth  decreases.  Conversely,  as  the 
bandwidth  of  the  jammer  approaches  infinity  the  time  domain  width  approaches  a  delta  function. 
The  null-to-null  width  of  the  time  domain  sa(a;)  (sine  over  argument)  curve  can  be  viewed  as  the 
time  in  which  the  statistics  of  the  jammer  remain  correlated.  Therefore,  if  the  time  to  transcend 
the  antenna  array  is  less  than  this  width  the  jammer  signal  remains  correlated.  In  other  words,  the 
relationship 


1 

Bi 


(3.37) 


must  hold  for  the  jammer  to  remain  correlated  across  the  array. 

The  jammer  signal  received  by  each  antenna  element  is  analyzed  in  the  same  manner  as  the 
target  signal.  Again,  the  ultimate  goal  is  development  of  a  space-time  covariance  matrix  R^*  for  use 
in  the  STAP  algorithm.  This  covariance  matrix  contains  the  correlation  statistics  of  the  jammer 
and  allows  the  algorithm  to  locate  the  jammer  spatially  and  in  time.  A  null  is  then  formed  in  that 
location. 
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Time  Domain  Representation 


Time 

Jammer  Power  Spectral  Density 


Frequency 

Figure  3.5  The  jammer  bandwidth  time  domain  correlation  relationship. 


First,  let  Xnm  represent  the  signal  received  from  the  jammer  on  the  antenna  element  at 
the  pulse.  Creating  a  vector  of  samples  across  the  antenna  array  elements  results  in 


Xm  — ■  OCfn^j  i&j )  )  j 


(3.38) 


where  the  jammer  amplitude  is  represented  by  am  and  is  the  jammer  spatial  steering  vector. 
Since  this  is  a  noise  jammer  the  amplitude  is  random.  This  amplitude  fluctuates  from  pulse-to-pulse 
and  can  be  represented  as  a  vector  of  random  amplitudes, 


aj  =  [ao  ai  ...  aM-if' ■ 


(3.39) 
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As  before,  the  Kronecker  Product  allows  the  compact  notation  of  a  space-time  snapshot  for  the 
jammer  Xj^ 

Xj  =  aj(8)aj.  (3.40) 


The  jammer  is  assumed  correlated  across  the  antenna  array  as  mentioned  in  the  beginning 
of  this  section.  Furthermore,  the  jammer  location  is  assumed  constant  throughout  a  CPI.  Using 
these  assumptions,  the  random  amplitudes  are  uncorrelated  as  shown  by 

where  6  is  again  the  Kronecker  delta  and  is  the  Jammer- to-Noise  Ratio  (JNR).  Writing  the  above 
equation  in  vector  form  results  in 

e{aa^}  =  (3.42) 


where  the  parameter  is  noise  power  per  element  as  mentioned  in  the  noise  model  section. 

The  JNR  is  easily  derived  using  the  radar  range  equation  as  given  by  Skolnik  [17].  First 
examine  the  one  way  range  equation  for  the  received  jammer  power 


Jo  = 


SjAe 


(3.43) 


given  Sj  (W/Hz)  is  the  radiated  power  spectral  density  (hence  it  already  includes  the  effect  of  the 
jammer  antenna  gain),  Ae  is  the  effective  area  of  the  receive  antenna,  and  Rj  is  the  range  to  the 
jammer.  The  effective  area  of  the  receive  antenna  depends  on  gain  and  wavelength  as  given  by 


47r  ’ 


(3.44) 
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where  Gr  is  the  gain  of  the  receive  antenna  and  \o  is  the  nominal  wavelength.  Substituting  for  the 
receive  antenna  effective  area  in  Eqn.  (3.43)  gives 


SjGAl 


(3.45) 


The  receive  antenna  elements  are  assumed  to  radiate  with  a  power  pattern  given  by  g{6y  or  the 
square  of  the  voltage  pattern  |/(0,  (^)p.  The  term  Lr  is  introduced  to  account  for  losses  in  the 
receiver  and  is  essentially  a  scaling  factor.  Using  these  tools,  the  received  jammer  power  spectral 
density  Jo  becomes 


Sj9{9,<P)Xl 

°  {4:T:fR]Lr 


(3.46) 


and  the  space-time  covariance  matrix  for  the  jammer  Rj  results  through  the  following  manipulation, 


=  £  {(oij  ®  &j){oti  (El 

=  •  (3-47) 

Again,  this  relationship  assumes  the  jammer  location  is  constant  over  the  CPI. 

Let  us  extend  Eqn.  (3.47)  to  two  jammer  signals.  The  vector  of  received  jamming  signals 
becomes 


Xrn  = 


(3.48) 
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where  the  superscripts  denote  the  jammer.  The  jammers  are  assumed  independent,  therefore  the 
cross  correlation  of  the  amplitudes  is  zero 

£  }  =  a^^f6i.v6ra-m'.  (3.49) 

Forming  the  space-time  snapshot  is  merely  a  matter  of  summing  the  signals  received  from  the  two 
jammers, 


(1) 

Xj  =  OL) 


I  -I-  af'*  I 


.(2) 


The  space-time  covariance  matrix  then  follows  directly, 


(3.50) 


R,  = 

=  R$‘>  +  Rf> 

=  I„  ®  (#5”  +  *f)  ,  (3.51) 

where  the  two  N  x  N  matrices  and  represent  the  jammer  spatial  covariance  matrices. 
They  are  the  result  of  taking  the  expected  value  of  the  received  jammer  signal  across  the  elements 
at  a  particular  pulse, 

^  ^  ~  •  (3.52) 

Now  take  a  look  at  the  more  general  case  of  P  jammers.  Since  the  result  is  again  a  sum  of 
the  individual  jammer  space-time  covariance  matrices,  Eqn.  (3.51)  still  holds  except  there  are  now 
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P  terms.  However,  it  can  be  more  compactly  represented  by  redefining  $  to  represent  all  of  them, 


$ 


A-SA^ 


(3.53) 


where  Aj  is  an  iV  x  iV  matrix  whose  columns  contain  each  of  the  jammers’  spatial  steering  vectors 
as  given  by 


=  K 


(1)  -(2) 


■n 


(3.54) 


and  Sj  is  a  P  X  P  matrix  containing  the  received  signal-to-noise  ratio  for  each  jammer  multiplied 
by  the  element  noise  power  along  its  diagonal,  defined  by 


3j  =  cr2diag(^j). 

Therefore,  the  jammer  space-time  covariance  matrix  for  P  jammers  becomes 


(3.55) 


Rj  =  Im  (3.56) 

3 .3 ,2,1  Jammer  Space-Time  Covariance  Matrix  Rank,  The  rank  of  the  jammer 
space-time  covariance  matrix  is  important  in  the  design  of  STAP  techniques.  If  the  matrix  is  of 
low  rank  the  STAP  algorithm  expends  fewer  degrees  of  freedom  placing  nulls  on  the  jammers. 

The  rank  of  the  space-time  covariance  matrix  for  a  single  jammer  is  equal  to  M  since  the 
outer  product  shown  in  Eqn.  (3.47)  is  a  linear  combination  of  the  vectors.  However,  for  P  jammers 
the  rank  increases  considerably. 

Since  the  size  of  Hj  is  MN  xMN  regardless  of  the  number  of  jammers,  the  worst  case  scenario 
is  a  matrix  of  full  rank  equal  to  {MNf,  However,  because  Rj  is  block  diagonal  (see  Eqn.  (3.56)), 
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the  rank  of  Rj  is  a  linear  combination  of  the  rank  of  each  of  the  blocks  along  the  diagonal  or 

rank(Rj)  =  M  rank($).  (3.57) 

The  problem  of  determining  the  jammer  rank  now  relies  on  the  rank  of  the  jammer  spatial  covariance 
matrix  The  following  relationship  between  the  dimension  of  two  matrices  holds  [19,  page  200], 

dim(A  4“  B)  =  dim(A)  +  dim(B)  —  dim(A  D  B).  (3.58) 

The  spatial  covariance  matrix  ^  representing  P  jammers  is  a  sum  of  the  individual  spatial 
covariance  matrices  (each  of  rank  one).  Given  statistically  independent  jammers  (different  loca¬ 
tions),  the  dimension  of  the  intersection  of  the  two  matrices  is  zero.  Therefore  the  dimension  (or 
rank)  of  the  jammer  spatial  covariance  matrix  ^  becomes 

rank(#)  =  P,  (3.59) 

where  P  is  the  number  of  jammers.  Using  this  fact,  the  rank  of  the  jammer  space-time  covariance 
matrix  is  seen  to  be  significantly  lower  than  full  rank, 

rank(Rj)  =  MP.  (3.60) 

Therefore  MP  Degrees  Of  Freedom  (DOF)  are  used  to  place  nulls  on  the  jammers  where  P  <  N. 

3,3.3  Clutter  Model  Clutter  is  defined  as  any  unwanted  radar  echo.  Because  the  focus 
of  this  section  is  on  airborne  surveillance  radars,  the  Earth  serves  as  the  primary  source  of  clutter. 
Hence,  the  model  described  in  this  section  only  describes  clutter  due  to  the  surface  of  the  Earth. 
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(a)  (b) 

Figure  3.6  Top  and  side  views  of  the  clutter  ring. 


The  return  from  clutter  is  very  similar  to  that  of  a  target.  However,  several  key  differences 
play  an  important  role  in  the  derivation  of  the  model.  First,  since  the  clutter  is  due  to  the  surface 
of  the  Earth  it  is  distributed  in  range  to  the  horizon.  Second,  the  surface  of  the  Earth  is  considered 
stationary.  Hence,  the  only  Doppler  shift  is  due  to  the  velocity  of  the  airborne  platform  relative  to 
the  clutter  patch  of  interest.  Using  these  observations,  the  clutter  is  broken  into  rings  of  constant 
range  around  the  aircraft.  Further  separation  of  the  clutter  ring  results  in  clutter  patches.  Each 
clutter  patch  is  considered  of  constant  range  and  velocity  with  respect  to  the  airborne  platform. 
Figure  3.1b  shows  the  clutter  ring  for  the  airborne  platform.  A  top  view  of  the  clutter  ring  (with 
the  airborne  platform  at  the  center)  is  shown  in  Fig.  3.6a. 

The  unambiguous  range  is  divided  into  L  range  bins  (rings  of  constant  range)  and  each  ring 
into  Nc  patches.  The  elevation  angle  to  the  clutter  ring/range  bin  is  given  by 


sin  ^ 


/ "I"  ha{ha  +  2ae)\ 

\  2Rc{(le  H"  ha)  ) 


(3.61) 


where  Rc  denotes  the  range  to  the  clutter  ring,  ha  the  altitude,  and  ae  the  effective  radius  of  the 
Earth.  As  is  customary,  the  effective  radius  is  considered  as  |  the  mean  radius  (|a)  of  the  Earth  to 
account  for  the  refraction  of  radio  waves  in  the  atmosphere  [17,  page  449].  As  defined  in  Fig.  3.1b, 
the  elevation  angle  is  referenced  to  the  horizontal  from  the  airborne  platform. 
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The  grazing  angle  is  the  angle  between  a  line  connecting  the  airborne  platform  to  the 
horizon  and  a  line  tangent  to  the  earth  at  the  horizon  point.  Shown  in  Fig.  3.6b,  this  angle  is  given 
by 


<0^  =  -  sin' 


/  —  ha{ha  +  2ae) 

V  2Rcae 


y 


(3.62) 


The  grazing  angle  is  zero  when  the  clutter  ring  of  interest  is  equal  to  the  range  to  the  radar  horizon. 

The  spatial  frequency  of  the  clutter  patch®  is  found  by  substituting  the  angles  associated 
with  the  patch  into  Eqn.  (3.15), 


= 


—  cos9ism<j>k, 
^0 


(3.63) 


where  Xq  is  the  nominal  wavelength,  9i  is  the  elevation  angle  to  the  patch,  and  (f>k  is  the  azimuth 
angle  to  the  patch.  The  normalized  Doppler  frequency  for  the  clutter  patch  is  expressed  as 

Since  the  clutter  is  extended  along  the  entire  Earth’s  surface  to  the  horizon,  the  unambiguous 
range  of  the  radar  becomes  an  important  consideration  in  the  development  of  the  model.  The 
radar’s  unambiguous  range  Ru  is  given  by  Skolnik  [17,  page  3]  as 


R 


u 


(3.64) 


where  c  is  the  propagation  velocity  and  fr  is  the  PRF.  Assuming  the  radar’s  unambiguous  range 
is  less  than  the  range  to  the  radar  horizon,  the  clutter  return  will  consist  of  echoes  from  multiple 

®The  clutter  patch  is  defined  as  the  patch  on  the  ring. 
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ranges.  Therefore,  each  ambiguous  range  is  given  by 


=  i?c  +  (^  -  ^)Ru  for  i  =  1, . . .  ,  Nr. 


(3.65) 


The  number  of  range  ambiguities  to  the  horizon  is  denoted  Nr. 

Prom  these  parameters,  the  component  of  the  space-time  snapshot  due  to  clutter  Xc  given 

by 

Nr  Nc  Nr  Nc 

Xc  =  Y^Y^O!ikv{’&ik,u>ik)  =  ^ ^ Q!i*b(a>ifc)  (gi  (3.66) 

i=l  fc=l  i=l  fc=l 

where  aik  denotes  the  random  amplitude  of  the  clutter  patch. 

The  next  step  is  to  determine  the  Clutter-to-Noise  Ratio  (CNR)  for  each  patch.  This  step 
is  similar  to  the  derivation  of  the  JNR.  However,  the  Radar  Cross  Section  (RCS)  of  each  clutter 
patch  must  be  determined  because  of  the  two  way  radar  equation.  The  dependence  of  the  RCS  of 
the  clutter  patch  is  explicit  to  highlight  the  fact  the  aircraft  is  flying  over  a  non-homogeneous 
surface  and  hence  the  RCS  per  unit  area  is  dependent  on  the  clutter  patch  under  consideration. 
The  RCS  per  unit  area  is  given  by  Skolnik  [17,  page  471]  as 


^k) 


Ac  ’ 


(3.67) 


where  ac  is  the  RCS  due  to  the  area  circumscribed  by  Ac.  Using  this  relationship,  the  RCS  of  the 
ik^^  clutter  patch  is  given  by 


^ik  ^  ^o(^ij0Aj)  X  Aik, 


(3.68) 
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where  the  Aik  is  the  area  of  the  clutter  patch  (assumed  constant  for  all  i  and  A;).  Using  the 
patch  area  given  by  Skolnik  [17,  page  472],  the  RCS  becomes 

cT 

o-ik  =  A(^-^sec'0i 

(* 

=  A0  2^  sec^i 

=  (^o{0iy(t>k)Ri  AR  secV^i.  (3.69) 


The  data  for  land  clutter  is  sometimes  given  by  the  parameter  7,  related  to  the  nominal  RCS 
of  the  clutter  patch  [17,  page  490]  through  the  relationship 


ao  =  ^sinipc,  (3.70) 


where  'ipc  is  again  the  grazing  angle.  This  model  serves  as  the  basis  for  the  clutter  simulations  found 
in  this  chapter.  The  terrain  dependent  parameter  is  taken  to  be  7,  however  Skolnik  states  7  also 
depends  on  the  grazing  angle  if  the  angle  of  incidence  is  near  perpendicular  or  the  grazing  angle  is 
extremely  low  in  the  case  of  a  low  flying  aircraft  [17,  page  490]. 

Using  the  above  deflnition  for  the  RCS  of  each  clutter  patch,  the  CNR  at  the  clutter 
patch  is  given  by  the  two  way  radar  range  equation. 


^ik  ^ 


PjTpGtiOij  (l>k^Q{0iy  (l>k)^^ik 

{Ai:fNoLsRi 


(3.71) 


where  the  element  pattern  is  represented  by  9{9i^<l>k)  as  before  and  the  two-sided  receiver  noise 
power  spectral  density  is  represented  by  such  that  the  noise  power  per  element  is  = 

^ofTp  (since  the  radar  is  assumed  to  have  a  separate  receiver  for  each  antenna  element). 
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Due  to  the  variability  of  the  clutter,  the  returns  from  different  patches  are  assumed  uncorre¬ 


lated  as  shown  by 


^  (3.72) 

Notice  the  element  noise  power  multiplied  by  the  CNR  results  in  the  power  in  the  clutter  patch 
as  expected.  Cross  correlation  terms  are  zero.  Using  this  fact,  the  clutter  space-time  covariance 
matrix  Rc  is  derived  as  follows. 


R.  = 


£{XcXf} 

{Nr  Nc  Nr  Nc 

i=lk=l  i'=lfc'=l  J 

{Nr  Nc  Nr  Nc  \ 

EEEE  aik  at'  fe-  ,  Wife)  ,Qi>k')^\ 

<=1  k=l  »'=1  *'=1  J 

Nt  Nc 


i=l  *=1 


(3.73) 


where  the  final  step  uses  the  results  from  Eqn.  (3.72).  The  clutter  space-time  covariance  matrix 
can  also  be  expressed  using  the  Kronecker  Product, 

Nr  Nc 

Rc  =  X ®  a(t?a)a(i?jfc)^.  (3.74) 

»=1  k=l 

An  even  more  compact  notation  is  available  through  the  definition  of  the  following  matrices,  similar 
to  that  done  for  the  jammer  covariance  matrix.  Begin  by  defining  the  steering  vector  matrix  with 
the  assumption  that  no  range  ambiguities  exist  {Nr  =  1), 


Vc  =  [vi;v2;--  -  ;vjvJ. 


(3.75) 


Then  form  the  clutter  power  distribution  into  the  diagonal  matrix 

S,  =  (72diag([6,6,---,^JvJ).  (3.76) 

Use  of  these  two  definitions  results  in  the  following  form  for  the  clutter  space-time  covariance 
matrix, 

Re  =  VeSeV".  (3.77) 


3,3,3. 1  Clutter  Ridges,  The  relationship  between  the  clutter  patch  and  the  motion 
of  the  airborne  platform  gives  rise  to  a  dependence  on  Doppler  frequency.  This  relationship  produces 
what  is  characterized  by  Ward  [21,  page  24]  as  a  clutter  ridge.  To  see  this,  first  revisit  the  spatial 
frequency  of  the  clutter  patch. 


'dik  = 

Ao 


(3.78) 


The  Doppler  frequency  of  this  patch  is  given  through  the  relationship  defined  previously  as 


fik  — 


2Vik 

Ao 


2k{ei,4>k)-Va 

Ao 


(3.79) 


where  Vik  is  the  velocity  of  the  clutter  patch  relative  to  the  airborne  platform  and  Va  is  the 
aircraft  velocity  vector.  If  the  aircraft  is  assumed  traveling  in  the  x  direction,  the  Doppler  frequency 
becomes 


fik  — 


2Va 

Ao 


cos^isin^fc. 


(3.80) 
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Rearranging  the  above  equation  highlights  the  expression  for  spatial  frequency, 


fik 


2dVa 

dXo 

d 


cos  6i  sin  (t>k 


(3.81) 


where  d  is  the  interelement  spacing  of  the  antenna.  Normalizing  this  expression  to  the  PRF  as 
before  results  in 


aik  =  ^=  fikTr  =  ,  (3.82) 

fr  ® 

where  Tr  is  the  PRI.  This  is  the  relationship  between  the  normalized  Doppler  frequency  and  the 
velocity  of  the  aircraft.  Since  the  spatial  frequency  of  the  clutter  patch  is  dependent  only  on 
(j)k  (the  aircraft  is  assumed  at  a  constant  altitude),  the  normalized  Doppler  frequency  is  linearly 
dependent  on  the  azimuth  with  a  slope  ^  given  by 


2VaTr 

d  ' 


(3.83) 


Examining  the  units  in  this  relationship  shows  it  represents  the  number  of  half-interelement 
spacings  d  crossed  during  one  PRI.  If  this  slope  is  greater  than  one  the  clutter  is  Doppler  ambiguous 
since  it  is  spread  over  a  frequency  band  wider  than  the  PRF.  If  it  is  less  than  one  it  is  unambiguous. 
To  further  understand  this  Doppler  ambiguity  relationship,  examine  the  Fourier  series  of  the  receive 
clutter  echo  pulse. 

The  Fourier  series  representation  of  a  rectangular  pulse  is  a  series  of  delta  functions  spaced  a 
distance  fr  apart  (the  PRF),  the  envelope  defined  by  a  sine  over  argument  curve.  The  Doppler  shift 
of  the  clutter  patch  return  moves  the  entire  spectrum  either  up  or  down  in  frequency,  depending 
on  whether  the  aircraft  is  approaching  or  receding  from  the  patch.  The  receiver  bandwidth  B  is 
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fixed  about  the  carrier  frequency  and  does  not  move  in  correlation  with  this  Doppler  shift.  The 
bandwidth  of  the  receiver  is  less  than  the  PRF.  Therefore,  if  the  Doppler  shift  is  enough  to  move 
the  spectrum  so  that  one  of  the  outlying  pulses  is  now  in  the  receiver  bandwidth,  the  radar  cannot 
distinguish  between  it  and  the  original  PRF  line.  Hence,  the  Doppler  frequency  is  ambiguous  or 
aliased. 

As  P  increases  it  becomes  more  difficult  to  suppress  the  clutter.  The  aliasing  in  the  Doppler 
frequency  space  explains  this  relationship.  Increased  aliasing  means  more  clutter  returns  contain 
the  same  Doppler  frequency  as  the  target,  further  obscuring  it.  Figure  3.18  shows  a  plot  of  an 
adapted  antenna  pattern  where  clutter  ridge  null  is  easily  visible  for  /3  =  1.  From  the  definition 
of  /?,  it  is  obvious  the  velocity  of  the  aircraft  plays  a  role  in  clutter  suppression.  Figure  3.7  shows 
several  plots  illustrating  the  aliasing  of  the  normalized  Doppler  frequency  for  various  values  of  /3. 

In  the  ^  =  1  plot  only  one  Doppler  frequency  exists  for  every  possible  spatial  frequency. 
Obviously,  for  values  of  P  less  than  one  this  relationship  still  exists.  This  is  termed  the  unambiguous 
Doppler  condition.  As  P  is  increased  more  and  more  spatial  frequencies  can  be  assigned  to  any  one 
Doppler  frequency.  Hence,  the  ambiguous  Doppler  condition.  Figure  3.7c  shows  a  case  where  as 
many  as  four  spatial  frequencies  can  be  associated  with  any  one  Doppler  frequency. 

3 *3,8,2  Rank  of  the  Clutter  Space- Time  Covariance  Matrix,  Because  the  clutter 
ridge  is  linear,  it  suggests  the  rank  of  the  clutter  matrix  is  low.  This  is  supported  by  Brennan^s 
rule.  This  rule  states 

rank(Rc)  »  [AT  +  (M  -  1)^J,  (3.84) 

where  the  [  J  symbols  indicate  rounding  to  the  nearest  integer.  Brennan’s  rule  states  the  rank  of 
the  clutter  covariance  matrix  is  a  function  of  the  number  of  antenna  elements,  the  CPI  length,  and 
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(c) 


(d) 


Figure  3.7  Clutter  ridge  plots  illustrating  the  normalized  Doppler  aliasing  for 
various  values  of  p. 
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Figure  3.8  Correlated  clutter  sampling  (/J  equal  to  an  integer)  gives  low  rank. 

the  slope  of  the  clutter  ridge.  For  integer  values  of  0  this  is  an  accurate  estimate  of  the  rank  as 
Fig.  3.11  shows  for  the  upcoming  example. 

Brennan’s  rule  makes  a  lot  sense.  It  is  primarily  a  function  of  the  number  of  half¬ 
interelement  spacings  transcended  during  one  PRI.  If  ^  is  an  integer  the  antenna  elements  basically 
take  the  same  clutter  sample  as  the  aircraft  travels.  Since  these  samples  correspond  to  the  same 
spatial  position,  the  clutter  space-time  covariance  matrix  rank  is  low  due  to  the  strong  correlation 
between  samples  (they  are  not  linearly  independent).  Figure  3.8  illustrates  this  concept. 

3. 3. 3,3  CNR  Plot  and  Clutter  Covariance  Matrix  Rank  Examples,  Given  the  el¬ 
ement  voltage  pattern  described  below  and  the  scenario  listed  in  Table  3.1,  let  us  calculate  the 
resulting  Clutter-to-Noise  Ratio  (CNR)  pattern  for  the  ik^^  clutter  patch  and  also  the  rank  of  the 
clutter  space-time  covariance  matrix  for  various  values  of  (3,  The  antenna  array  is  rectangular  with 
elements  in  the  x  and  z  directions.  The  beamforming  operation  uses  only  a  row  of  elements  in  the 
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Table  3.1  Scenario  parameters  for  the  clutter  example. 


Variable 

Value 

M 

18 

N 

18 

fo 

450  MHz 

fr 

300  Hz 

Tp 

200  fjs 

Pi 

200  kW 

B 

4  MHz 

Fn  (Noise  Figure) 

3dB 

Nc 

361 

h/Q 

9  km 

Rc 

130  km 

7 

-3  dB 

Transmit  Gain 

22  dB 

Column  Receive  Gain 

10  dB 

Number  of  Elements  x 

18 

Number  of  Elements  z 

4 

Element  Pattern 

Cosine 

Element  Gain 

4dB 

Element  Backlobe  Level 

-30  dB 

Transmit  Taper 

Uniform  (None) 

X  direction,  18  elements,  each  with  the  voltage  pattern 


cos  (j)k  —90®  <(t>k<  90®  frontlobe 

< 

be  COS  <l>k  90®  <<t>k<  270®  backlobe. 


The  value  of  he  required  to  give  a  backlobe  of  -30  dB  is  0.032,  the  inverse  of  the  dB  operation. 
Figure  3.9  shows  the  antenna’s  pattern  using  this  value. 

The  transmit  gain  in  Eqn.  (3.71)  is  a  function  of  the  azimuth  and  elevation  angles  as  given 
by  the  array  pattern.  The  array  pattern  is  calculated  using  the  same  techniques  as  in  Chapter  II 
and  is  given  by 

i=0 
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Antenna  Element  Pattern  (dB) 


270 


Figure  3.9 


Element  voltage  pattern  for  the  clutter  example. 


where  k  is  the  spatial  frequency  vector  and  Wi  is  the  individual  antenna  element  weight.  No  weights 
have  been  applied  to  the  antenna  element  signals,  hence  w  is  equal  to  one  for  all  elements.  Also, 
the  antenna  beam  is  not  steered  in  this  example  so  ko  is  equal  to  the  zero  vector,  resulting  in  a 
beam  along  the  boresight  of  the  array.  Using  Eqn.  (2.3),  the  spatial  frequency  is  given  by 

0<yr  /\ 

k  = 

2‘jj’ 

=  —  {cos  6i  sin  (l>kX  +  cos  6i  cos  (f)ky-^  sin  6iz} . 

Ao 

Notice,  given  the  antenna  elements  oriented  in  the  x  direction,  the  formula  coincides  with  the 
expression  for  spatial  frequency  found  in  Eqn,  (3.15),  Since  the  beamforming  operation  is  only 
across  the  elements  in  the  x  direction,  the  outer  product  in  the  exponential  simplifies  the  equation 
to 


N-l 


n=0 


The  summation’s  argument  can  be  written  to  the  power  of  it’s  index  and  simplified, 


W{eu(j>k) 


Jl  ^  ^g-J^cos04Sin</)A, 

^  n=0 

sin  ( cos  6i  sin 
sin  <f>k  V  ^  / 

sin  cos  di  sin  j 


N 


The  transmit  gain  pattern  for  the  array  is  now  given  by  (the  antenna  pattern  as  derived  in  Chapter  II 
is  voltage) 


Gti9i,h)  =  9{ei,h)W\0i,cl>k), 
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where  the  element  power  pattern  9{0i,<l)k)  is  simply  the  voltage  pattern  squared,  \f{9i^(l>k)\^. 

The  next  step  is  to  determine  the  elevation  angle  for  the  given  scenario.  Using  Eqn.  (3.61) 
and  the  values  given  in  Table  3.1,  the  elevation  angle  is  calculated  as  -0.0769  radians.  The  cosine 
of  this  angle  is  approximately  one  and  is  ignored  from  this  point  forward. 

The  next  major  calculation  is  the  determination  of  the  RCS  from  the  clutter  patch.  To 
begin,  determine  the  grazing  angle  for  the  aircraft’s  given  altitude  and  range  of  interest.  Substi¬ 
tuting  the  appropriate  values  into  Eqn.  (3.62)  shows  the  grazing  angle  to  be  0.0617  radians.  Prom 
here,  the  RCS  per  unit  area  is  easily  calculated  given  the  terrain  parameter  7, 

(To  =  10“°*®sin(0.0617)  =  0.0309. 


The  granularity  of  the  patches  around  the  clutter  ring  of  interest  is  given  by 


^  =  0.0174  radians 

■Af (*  361 


and  the  radar’s  range  resolution  is 


.  „  c  3  X 10® 

A/?  =  —  =  — TT7T  =  37.5  m. 


2B  2(4  X  10®) 


This  example  is  performed  at  the  first  ambiguous  range,  therefore  =  Also,  our  interest 
is  focused  on  only  one  clutter  ring.  Using  these  assumptions  and  the  constant  7  clutter  model,  the 
RCS  of  the  ifc**  clutter  patch  is  constant  around  the  clutter  ring  and  equal  to 


<^ik  =  (ToRiA^Aiisec^j 


=  0.0309  X  (130  X  10®)  X  0.0174  x  37.5  x  sec(0.0617) 
=  2626.09  m®. 
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The  last  unknown  parameter  for  the  CNR  pattern  (Eqn.  (3.71))  is  the  receiver  noise  power 
spectral  density,  Ng.  The  system  noise  figure  is  given  in  the  table.  Therefore,  the  noise  power 
spectral  density  in  the  receiver  is  [17,  page  19] 

Fn  =  or  No  =  FnkToGr, 

where  the  gain  of  the  receiver  Gr  is  assumed  equal  to  one,  the  system  noise  figure  Fn  is  given  as 
part  of  the  problem,  k  is  Boltzman’s  constant,  and  To  is  the  standard  temperature  (290®K).  Solving 
for  No  gives 


No  =  10°'^  X  (13.80658  x  lO"^^)  x  290  =  319.554  x  10■^^ 

This  is  the  noise  due  to  thermal  activity  in  the  receiver. 

Substituting  these  values  and  an  interelement  spacing  d  equal  to  half  the  wavelength  into 
Eqn.  (3.71)  results  in  the  plot  shown  in  Fig.  3.10.  This  is  the  CNR  around  the  entire  clutter  ring. 
Since  the  antenna  beam  is  not  steered,  the  main  beam  is  located  at  zero  degrees  azimuth.  Notice 
the  sharp  roll-off  at  ±90°  azimuth  due  to  the  antenna  element  patterns.  These  locations  represent 
the  ends  of  the  antenna  array. 

The  next  part  of  the  example  examines  the  performance  of  Brennan’s  rule  in  determining  the 
rank  of  the  clutter  space-time  covariance  matrix  Rc.  Since  we  want  to  plot  the  rank  of  the  clutter 
matrix  for  various  values  of  first  determine  the  Doppler  frequency  of  the  clutter  patch  in 
terms  of  this  using  Eqn.  (3.83), 


fik  — 


^  cos  0i  sin 


2l3d 

2TrXo 

TrXn 


cos  0i  sin 


cos  9i  sin 
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Figure  3.10  CNR  around  the  clutter  ring  (main  beam  pointed  at  0°  azimuth). 
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The  temporal  steering  vector  requires  the  normalized  Doppler  of  the  clutter  patch  (see  Eqn.  (3,27)) 
given  by 


sin  <f>k- 


Prom  this  point  on,  the  formulation  of  the  clutter  matrix  is  the  implementation  of  Eqn.  (3.74)  in 
Matlab®.  Because  of  the  peculiarities  of  Matlab®  the  implementation  is  not  discussed  here. 
The  element  noise  power  is  assumed  equal  to  one  for  this  simulation. 

Figure  3.11  shows  the  resulting  ranks  for  clutter  matrices  corresponding  to  /3  values  of  0.6,  1, 
2,  2.83,  and  3.  Because  the  number  of  nonzero  eigenvalues  is  synonymous  with  the  rank  of  a  matrix, 
this  metric  is  used.  As  mentioned  previously,  Brennan’s  rule  predicts  the  rank  of  the  matrix  very 
well  when  ^  is  an  integer.  This  means  the  number  of  half-interelement  spaces  transcended  during 
a  PRI  is  an  integer  value.  Therefore  the  antenna  elements  line  up  and  the  subsequent  elements 
basically  retrieve  the  same  clutter  sample  a-s  the  first.  Hence,  the  clutter  matrix  contains  correlated 
clutter  samples.  This  condition  is  shown  in  Fig.  3.8.  However,  Brennan’s  rule  does  not  hold  when 
given  uncorrelated  clutter  samples  (non-integer  values  of  P)  as  shown  by  the  tapering  eigenvalues 
for  the  p  =  0.6  and  ^  =  2.83  curves. 


34  Direct  Form  Fully  Adaptive  STAP 

The  purpose  of  Space-Time  Adaptive  processing  is  to  direct  an  antenna  null  in  the  direction 
of  a  jammer  or  a  significant  clutter  source.  This  minimizes  the  effect  of  the  interference  and  can 
lead  to  the  detection  of  a  target  in  environments  standard  radar  processing  cannot  handle.  This 
section  discusses  computation  of  the  antenna  pattern/filter  weights  to  steer  the  antenna  nulls  using 
a  fully  adaptive  algorithm. 

The  two  basic  approaches  to  STAP  are  Fully  Adaptive  STAP  and  Partially  Adaptive,  or 
Reduced-Dimension  STAP.  There  are  inherent  advantages  and  disadvantages  to  both.  Obviously, 
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Eigenvalue  Number 


Figure  3.11  Illustration  of  Brennan’s  rule  for  clutter  matrix  rank  using  various 
values  of  /?. 


the  best  performance  to  be  obtained  is  with  a  fully  adaptive  algorithm.  However,  the  computational 
burden  becomes  too  much  to  perform  within  the  timeline  of  the  typical  radar. 

Fully  adaptive  STAP  results  in  a  weight  vector  for  every  element  of  the  space-time  snapshot 
(see  Eqn.  (3.31)),  i.e.  every  element  and  pulse  in  the  CPI  is  weighted  by  the  digital  processor.  This 
entire  process  is  described  by 

x  =  w^x-  (3.85) 

The  vector  w  is  the  weight  vector  defined  by  [21] 

w  =  R-^vt,  (3.86) 

where  R  is  defined  as  the  sum  of  all  the  undesired  component  space-time  covariance  matrices, 
Yi  is  the  target  steering  vector,  and  z  is  the  output  of  the  algorithm.  The  undesired  space-time 
covariance  matrix  is  given  by 


R  ~  Ry^  -j-  Rj  +  Rcj 


(3.87) 


where  R„  represents  thermal  noise,  R^  represents  jammers,  and  Rc  represents  clutter. 

The  expression  for  the  algorithm  (Eqn.  (3.85))  maximizes  Signal-to-Interference-plus-Noise- 
Ratio  (SINK),  maximizes  probability  of  detection  for  a  given  false  alarm  probability,  and  with  the 
proper  choice  of  scale  factor  minimizes  the  output  noise  power  subject  to  a  unity  gain  constraint  in 
the  target  direction  [21,  page  57].  One  problem  with  this  approach  is  high  sidelobes  in  both  angle 
and  Doppler.  For  this  reason  Ward  also  introduces  the  tapered  weight  vector, 


gt  =  tOVf, 


(3.88) 
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where  t  is  a  two  dimensional  window/taper  and  0  is  the  Hadamard  product  (equivalent  to  an 
element  by  element  multiply).  This  two  dimensional  window  can  be  generated  from  two  one  di¬ 
mensional  windows  using  the  Kronecker  product  as  shown  by 

t  =  tft(8)ta,  (3.89) 

where  tj,  and  to  are  the  one  dimensional  window  functions.  An  expression  equivalent  to  Eqn.  (3.88) 
is 


gt  =  (t6  0b<)0(to0at). 


(3.90) 


This  tapered  space-time  steering  vector  can  be  substituted  for  v<  in  Eqn.  (3.86).  The  window 
in  this  new  weight  vector  produces  low  sidelobes  at  the  sacrifice  of  mainbeam  width  and  gain. 
Correspondingly,  it  no  longer  produces  the  optimum  SINK.  Harris  [11]  offers  an  excellent  discussion 
on  windowing  with  a  direct  comparison  of  the  different  types  of  windows.  The  SINK  is  defined 
as  [21] 


SINR  = 


w^Rw  ’ 


(3.91) 


where  the  numerator  represents  the  signal  power  and  the  denominator  represents  the  interfer¬ 
ence  plus  noise  power.  To  get  the  optimum  SINK,  substitute  the  optimum  weight  vector  w  from 
Eqn.  (3.86)  into  the  above  expression  and  simplify, 


SINRopt 


vfR-iRR-ivt 

(T^gtvfR-^VtvfR-^Vt 

vfR-ivi 


(3.92) 
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where  cr^  is  the  noise  power  per  element  and  is  the  target  SNR  per  element.  Since  Eqn.  (3.92)  is 
for  one  particular  target  azimuth  and  Doppler,  it  can  be  expanded  across  normalized  Doppler  by 
highlighting  the  temporal  steering  vector  dependence  on  a), 


SINRopt(ci^)  =  (bt(c4;)  (8)  R  ^  (bt(c<;)  <8)  at) . 


(3.93) 


The  equation  could  also  be  written  to  highlight  the  dependence  on  spatial  frequency  because  of  the 
spatial  steering  vector  but  then  the  transmit  and  receive  antenna  patterns  have  to  be  taken  into 
account  when  plotting  this  expression. 

Similarly,  the  SINK  for  the  tapered  pattern  can  be  obtained  by  substituting  Eqn.  (3.88)  into 
Eqn.  (3.91)  and  simplifying. 


SINRtap  = 


cr^^t  |gf 

gfR-ig< 


(3.94) 


Highlighting  the  dependence  on  normalized  Doppler  a)  produces 


SINRtap(a>)  = 


2 

[(*6  ©  bt(a>))  ®  (to  ©  ai)]'^  R”’-  (bf (w)  0  a.t) 

[(tj  ©  bt(<I;))  0  (to  ©  at)]'*^  R-^  [(tj  0  bt(w))  0  (to  ©  a*)] 


(3.95) 


Figure  3.12  shows  the  effect  of  three  different  types  of  windows  on  the  SINK  for  a  fully 
adaptive  processor.  The  plots  shown  are  for  an  18  azimuth  element  array  with  18  pulses  in  the 
CPI.  The  noise  power  per  element  and  target  SNR  per  element  per  pulse  are  both  set  to 
one  for  the  plots.  The  interference  matrix  R  is  the  same  as  that  generated  for  the  previous  clutter 
example  with  the  addition  of  two  jammers  and  noise.  The  windows/tapers  were  applied  in  both  the 
spatial  and  temporal  domains.  Referencing  Pig.  3.12  shows  the  poorest  performance  of  the  three 
tapers  shown  is  with  the  Blackman  Harris  window.  This  poor  performance  is  expected  because 
of  the  extremely  low  -90  dB  sidelobes  produced  by  this  window  in  both  the  spatial  and  temporal 
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Fully  Adaptive  STAP  Tapered  Steering  Vector  SINR  Comparison 


Figure  3.12  An  SINR  comparison  of  fully  adaptive  STAP  for  the  Blackman  Har¬ 
ris,  Hanning,  and  Hamming  [11]  tapered  steering  vectors. 


domains.  The  loss  in  main  beam  gain  outweighs  the  gain  found  by  suppressing  the  interference 
with  the  extremely  low  sidelobes. 

The  SINR  represents  a  good  tool  for  the  comparison  of  various  STAP  algorithms.  Another 
tool  is  the  SINR  loss.  This  metric  compares  the  performance  of  the  algorithm  to  the  optimal  case,  a 
noise  only  environment.  In  the  noise  only  environment,  the  interference  covariance  matrix  becomes 


R  =  Rn  = 


(3.96) 
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Substituting  this  noise  only  covariance  matrix  into  Eqn.  (3.92)  results  in  the  optimum  Signal-to- 
Noise  Ratio  (SNR), 


SNRopt  =  ^Imjvv<  =  6vf  vt.  (3.97) 

Since  the  space-time  steering  vector  is  the  complex  exponential  form  of  Eqns.  (3.27)  and  (3.26), 
it’s  inner  product  is  simply  the  number  of  elements  in  the  vector,  or  MJV.  Therefore,  the  optimum 
SNR  becomes 


SNRopt  =  MN^t, 


(3.98) 


where  is  the  target  SNR  per  element  per  pulse.  Using  this  optimum  SNR,  we  can  define  the 
SINR  loss  of  the  space-time  processing  algorithm  as 


I'SINR  = 


SINR(^) 

SNRopt  ■ 


(3.99) 


A  plot  of  this  metric  corresponding  to  the  windows  used  in  Fig.  3.12  is  shown  in  Fig.  3.13. 

A  block  diagram  of  the  basic  STAP  process  is  shown  in  Fig.  3.14.  The  process  involves  three 
major  steps: 


•  Train  the  processor,  i.e.  generate  a  second  set  of  data  for  use  in  the  computation  of  the  filter 
weights.  The  entire  CPI  datacube  is  used  for  the  training  algorithm  in  fully  adaptive  STAP. 

•  Generate  the  weights  using  the  new  set  of  data  obtained  from  step  one. 

•  Apply  the  weights  to  the  space-time  snapshot  at  the  range  gate  of  interest. 

Once  the  weights  have  been  applied,  the  result  is  sent  to  the  radar  threshold  detection  circuit 
to  make  a  target  decision.  A  more  detailed  diagram  of  the  fully  adaptive  STAP  process  is  shown  in 
Fig.  3.15.  The  enormous  amount  of  data  calculation  is  more  readily  apparent  from  this  diagram. 
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Fully  Adaptive  STAP  Tapered  Steering  Vector  SINK  Loss  Comparison 


Figure  3.13  An  SINK  loss  comparison  of  fully  adaptive  STAP  for  the  Blackman 
Harris,  Hanning,  and  Hamming  [11]  tapered  steering  vectors. 
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Figure  3.14  Top  level  block  diagram  for  a  space-time  adaptive  processor. 
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Figure  3.15  Fully  adaptive  space-time  processing  model. 


Another  metric  of  a  STAP  algorithm’s  performance  is  the  antenna  pattern  it  produces.  The 
adapted  pattern  is  given  by  [21] 

Pw('i?,w)  =  |w^v(i9,fD)p,  (3.100) 

where  w  is  the  weight  vector  and  v  is  the  steering  vector  (not  necessarily  steered  in  the  direction 
of  the  target). 

An  example  of  this  adapted  pattern  using  the  parameters  from  the  clutter  example  is  shown 
in  Figs.  3.16,  3.17,  and  3.18.  These  patterns  also  show  the  introduction  of  two  jammers,  both  at 
the  same  elevation  as  the  airborne  platform.  The  azimuth  locations  of  the  jammers  are  —40®  and 
25®  at  a  range  of  370  km.  The  effective  radiated  power  density  for  each  jammer  is  1  kW/MHz.  The 
target  is  located  at  0®  azimuth  with  a  Doppler  frequency  of  100  Hz.  The  value  of  /?  for  the  plots  is 
1,  reflected  by  the  clutter  ridge  null  shown  in  the  fully  adapted  pattern  found  in  Fig.  3.18. 
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Adapted  Pattern 


Figure  3.16  Fully  adapted  antenna  pattern  at  lOOHz  Doppler  frequency.  Vertical 
lines  indicate  jammer  locations. 
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Adapted  Pattern 


Figure  3.16  shows  a  cut  along  the  target  Doppler  from  Fig.  3.18.  Since  the  target  is  located 
at  0^  azimuth,  the  mainbeam  of  the  pattern  is  located  at  0®  azimuth.  Notice  the  nulls  placed 
along  the  jammer  locations  (sin(— 40°)  =  —0.643  and  sin(25°)  =  0.423)).  The  nulls  actually  extend 
much  lower.  The  plot  doesn’t  show  this  because  of  the  granularity  of  the  azimuth  angles  used  to 
generate  it.  Calculations  at  the  these  exact  azimuth  values  show  the  nulls  possess  greater  than  90 
dB  attenuation. 

A  cut  along  the  target’s  azimuth  is  shown  in  Fig.  3.17.  The  mainbeam  is  placed  at  the  target’s 
Doppler  at  100  Hz  with  a  strong  null  at  zero  Doppler.  Zero  Doppler  corresponds  to  mainbeam 
clutter,  a  significant  source  of  target  acquisition  difficulties.  Again,  the  granularity  of  the  Doppler 
frequency  values  used  to  generate  the  plot  does  not  show  the  true  depth  of  the  null.  A  calculation 
at  exactly  zero  Doppler  shows  this  null  is  greater  than  90  dB  down  as  in  the  case  of  the  jammers. 
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Figure  3.18  Three  dimensional  image  of  the  fully  adapted  pattern  Py,  with  no 
tapering  for  a  single  target  and  two  jammers. 


The  final  plot  (Fig.  3.18)  shows  the  entire  pattern  as  viewed  from  the  surface.  The  color 
represents  the  amplitude  of  the  lobes.  This  figure  shows  the  clutter  ridge,  corresponding  to  /?  =  1. 
Notice  the  ridge  corresponds  to  the  one  predicted  by  Fig.  3.7c.  It  is  not  linear  in  this  plot  because 
it  is  plotted  against  ^  instead  of  sin  <^. 

Fully  adaptive  STAP  is  too  computationally  intensive  to  implement  real  time.  The  number 
of  computations  is  on  the  order  of  (MN)^  [21].  Real  time  calculation  is  beyond  current  technology, 
although  the  computer  industry  is  improving  processor  speed  in  such  large  steps  that  the  ability  to 
implement  this  algorithm  could  be  reached  in  the  near  future  (possibly  within  the  next  few  years). 
Fully  adaptive  STAP  is  used  primarily  as  the  measure  to  compare  reduced  dimension  algorithms. 
Reduced  dimension  algorithms  can  be  implemented  within  the  constraints  of  current  computational 
technology. 
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Figure  3.19  The  Factored  Time  Space  (FTS)  reduced  dimension  algorithm  block 
diagram.  This  diagram  applies  to  one  range  gate  and  is  repeated  for  every  range 
gate  of  interest. 


5.5  Partially  Adaptive  Factored  Time-Space  STAP 

As  just  discussed,  computational  load  drives  the  need  for  a  reduced  dimension  algorithm. 
One  of  the  simplest  methods  is  the  Factored  Time-Space  (FTS)  algorithm. 

This  algorithm  reduces  the  computational  load  by  sacrificing  the  Doppler  Degrees  Of  Freedom 
(DOF).  The  incoming  signal  from  each  antenna  element  is  first  Doppler  filtered  using  the  pulse 
returns  from  the  range  gate/cell  of  interest.  A  noise  covariance  matrix  is  then  estimated  from  the 
Doppler  filtered  data  and  used  to  adaptively  filter  in  the  spatial  domain,  resulting  in  a  scalar  output 
for  each  Doppler  bin.  A  block  diagram  of  this  process  is  shown  in  Fig.  3.19.  The  operations  shown 
in  this  block  diagram  are  repeated  for  each  range  cell  of  interest. 

Although  Fig.  3.19  shows  a  Doppler  filter  bank  of  M  bins,  equal  to  the  number  of  pulses  in 
the  CPI,  the  number  of  bins  does  not  in  general  have  to  be  equal  to  the  number  of  pulses.  The 
Doppler  filter  is  implemented  with  the  Discrete  Fourier  Transform  (DFT).  The  use  of  the  DFT 
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allows  us  to  zero  pad  the  incoming  data  stream  (pulse  returns  from  a  particular  antenna  element) 
to  obtain  a  greater  number  of  Doppler  bins. 

Using  the  DFT  matrix  allows  us  to  mathematically  express  the  FTS  algorithm.  Since  this 
algorithm  results  in  a  single  scalar  for  each  Doppler  bin,  let’s  consider  the  mathematical  architecture 
behind  just  a  single  Doppler  bin,  m.  Begin  by  first  obtaining  the  Doppler  filter  bank  output  for 
this  bin.  First  define  the  Doppler  filter  bank  using  the  M  x  Af  DFT  matrix  W.  Since  we  have 
sacrificed  all  the  DOF  in  the  Doppler  domain,  we  need  a  Doppler  filter  bank  with  very  low  sidelobes 
to  suppress  clutter.  Therefore,  let  us  express  this  low  sidelobe  Doppler  filter  bank  as 

F  =  diag(t6)W^,  (3.101) 


where  tj,  is  a  window  function  (see  Harris  [11]  for  a  complete  discussion  of  windowing)  and  diag 
is  the  Matlab®  function  for  creating  a  diagonal  matrix  from  the  vector  argument.  The  diagonal 
elements  of  the  matrix  are  the  elements  of  tj.  Breaking  this  windowed  DFT  matrix  down  into  its 
columns, 


F  = 


fo  fi 


fM-l  » 


(3.102) 


allows  us  to  write  the  Doppler  filter  bank  output  for  a  particular  bin.  This  output  for  Doppler  bin 
m  becomes 


Xm  =  (fm®Ijv)"x,  (3.103) 

where  x  is  the  space- time  snapshot  from  Eqn.  (3.29).  The  spatially  adaptive  filter  weights  for  this 
Doppler  bin  become 

Wm  =  Rm  (3.104) 
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where  a*  is  the  spatial  steering  vector  as  defined  in  Eqn.  (3.26)  (or  a  tapered  spatial  steering  vector 
gt  could  be  substituted)  and  the  covariance  matrix  is  given  by 


^  f  {XmXm}  =  (fm  ®  liv)"  R(fm  ®  Ijv)  •  (3.105) 

Using  these  weights,  the  output  of  the  FTS  algorithm  for  the  Doppler  bin  Zm  is 

Zm  =  W^Xm-  (3.106) 

If  we  substitute  Eqn.  (3.103)  into  Eqn.  (3.106),  an  expression  results  using  the  full  dimensioned 
incoming  data  vector  or  space-time  snapshot  Xj 

^m  =  w"(f^®I^)"x-  (3.107) 

A  full  dimension  weight  vector  associated  with  the  above  equation  can  be  defined  as 

w"  =  w"(f^®Ijv)"  (3.108) 

or 

=  (fm  ^  liv)  Wm  =  fm  0  Wm«  (3.109) 

This  gives  the  alternate  form  of  the  adaptive  processor’s  output, 

Zm  =  yfmX-  (3.110) 
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We  can  use  the  SINR  as  a  comparison  of  the  FTS  algorithm  to  fully  adaptive  STAR  Substi¬ 


tuting  Eqn.  (3.109)  into  Eqn.  (3.91)  results  in  an  SINR  for  Doppler  bin  m, 


SINRm  = 


(fm<S)w^)^(b<(a;)(8)a<) 


(fm  ®  WjTi)  R(fyj^0Wni) 


(3.111) 


This  equation  can  be  further  simplified  by  substituting  the  FTS  weight  vector  from  Eqn.  (3.104) 
into  Eqn.  (3.111), 


SINR„i  = 


(3.112) 


where  gt  is  the  tapered  spatial  steering  vector, 


gt  =  taOat. 


(3.113) 


Figure  3.20  compares  the  SINR  Loss  Lsinr  for  the  FTS  algorithm  using  several  different 
Doppler  filter  bank  tapers,  t^.  No  taper  has  been  used  in  the  spatial  domain.  Since  there  are  18 
pulses  in  the  CPI,  the  Doppler  filter  contains  18  bins.  There  is  an  SINR  expression  for  each  bin. 
The  figures  plot  the  envelope  of  all  the  SINR  expressions  for  each  bin.  As  the  figure  illustrates, 
the  lower  sidelobes  in  the  Doppler  filter  bank  increase  the  SINR.  This  trend  is  expected  since  the 
primary  benefit  of  the  lower  sidelobes  is  suppression  of  the  interference. 

Figure  3.20  also  brings  forth  an  unexpected  point.  As  the  quality  of  the  Doppler  filter  increases 
(as  measured  by  the  lower  sidelobes  due  to  the  window  function  used),  the  performance  of  the  FTS 
algorithm  approaches  the  optimum  for  fully  adaptive  STAP.  The  non-adaptive  Doppler  filter  is 
performing  almost  as  well  as  the  adaptive  Doppler  filter.  The  trade-off  occurs  on  the  edges  of  the 
normalized  Doppler  spectrum.  Considerable  roll-off  occurs  above  Q  =  0.8  and  below  Q  =  0.2  for 
the  Blackman  Harris  Doppler  window. 
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FTS  SINK  Loss  Comparison 


Figure  3.20  Comparison  of  Lsinr  for  the  FTS  algorithm  with  four  Doppler  Filter 
bank  windows  and  no  spatial  taper.  The  optimum  SINR  curve  corresponds  to  a 
fully  adaptive  algorithm  (18  x  18  =  324  DOF)  with  no  tapers.  The  FTS  algorithm 
is  adaptive  in  the  spatial  domain  only,  resulting  in  18  DOF  for  this  example.  The 
Blackman  Harris  window  results  in  the  greatest  SINR  for  the  FTS  algorithm  because 
of  the  -92  dB  sidelobes.  The  Hamming  and  Hanning  window  sidelobes  are  -43  dB 
and  -35  dB,  respectively.  The  greater  sidelobe  roll-off  of  the  Hanning  window  (-18 
dB/Octave)  allows  it  to  outperform  the  Hamming  window  (-6  dB/Octave).  Window 
parameters  taken  from  Harris  [11]. 
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FTS  SINR  Comparison 


Figure  3.21  Comparison  of  Lsinr  for  the  FTS  algorithm  with  four  Doppler  Filter 
bank  windows  and  the  same  spatial  window  as  the  Doppler  window.  The  introduc¬ 
tion  of  the  spatial  taper  reduced  the  SINK  as  expected.  The  optimum  SINR  curve 
corresponds  to  a  fully  adaptive  algorithm  (18  x  18  =  324  DOF)  with  no  tapers. 
The  FTS  algorithm  is  adaptive  in  the  spatial  domain  only,  resulting  in  18  DOF  for 
this  example.  The  Blackman  Harris  window  used  both  in  the  non-adaptive  Doppler 
filter  and  the  spatial  taper  results  in  the  greatest  SINR  for  the  FTS  algorithm  be¬ 
cause  of  the  -92  dB  sidelobes.  The  Hamming  and  Hanning  window  sidelobes  are  -43 
dB  and  -35  dB,  respectively.  The  greater  sidelobe  roll-off  of  the  Hanning  window 
(-18  dB/Octave)  allows  it  to  outperform  the  Hamming  window  (-6  dB/Octave). 
Window  parameters  taken  from  Harris  [11]. 


As  a  final  comparison,  a  spatial  filter  matching  the  Doppler  filter  was  added  to  the  FTS 
algorithm.  These  SINR  Loss  plots  are  shown  in  Fig.  3.21.  Each  curve  shows  a  decrease  in  SINR 
due  to  the  addition  of  the  spatial  taper. 


3.6  Summary 

This  chapter  introduced  the  basic  fundamentals  of  fully  adaptive  STAR  The  fully  adaptive 
algorithm  discussed  placed  a  null  at  zero  degrees  azimuth  to  suppress  mainbeam  clutter.  It  also 
placed  nulls  on  any  barrage  noise  jammers  present.  However,  one  point  not  addressed  is  how  the 


3-55 


radar  informed  the  operator  of  the  position  of  the  antenna  nulls.  If  the  algorithm  did  a  very  good 
job  of  filtering  out  the  noise  jammers,  it  seems  theoretically  possible  for  the  radar  operator  to 
receive  no  indication  of  any  targets  coming  in  at  that  azimuth/elevation  (they  would  be  filtered 
out).  An  obvious  ECCM  tactic  would  be  to  equip  an  aircraft  with  a  barrage  noise  jammer.  A  radar 
equipped  with  this  algorithm  would  never  know  the  target  was  there. 

Other  items  of  interest  were  the  series  of  nulls  in  the  antenna  pattern.  These  graphs  not 
only  showed  the  deep  antenna  nulls  at  jammer  azimuths  and  zero  Doppler  (mainbeam  clutter)  but 
also  a  series  of  nulls  (13  of  them)  across  the  azimuth  range.  These  nulls  range  from  -50  to  -70  dB 
compared  to  the  -90  of  the  jammer  nulls.  These  nulls  not  only  filter  out  undesired  returns  but  also 
any  desired  returns  at  these  azimuths.  The  Doppler  pattern  suffers  the  same  problem. 

Because  of  the  computational  issues  associated  with  the  fully  adaptive  STAP  algorithm, 
the  partially  adaptive  FTS  algorithm  was  introduced.  There 
algorithms  available,  this  one  was  introduced  because  of  its 
reduced  dimension  concepts. 


are  many  other  reduced  dimension 
simplicity  and  value  in  explaining 
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IV.  Reduced  Dimensionality  Through  the  Cross  Spectral  Metric  Algorithm 

4.1  Introduction 

As  mentioned  previously,  Ward’s  [21]  work  serves  as  an  excellent  tutorial  for  STAP.  Chap¬ 
ter  III  is  a  summary  of  this  work,  serving  to  introduce  the  basics  of  STAP  and  develop  an  under¬ 
standing  of  the  important  issues.  This  next  chapter  delves  deeper  into  the  knowledge  required  to 
meet  the  main  objective  of  the  thesis:  determination  of  the  secondary  data  support  required  for  one 
of  the  recent  approaches  to  reduced  dimension  STAR  A  secondary  objective  involves  determining 
the  effect  of  non-homogeneities  in  the  secondary  data  of  this  model.  This  objective  also  requires 
the  information  found  here. 

The  Cross  Spectral  Metric  (CSM)  algorithm  just  recently  introduced  by  Goldstein  [6]  was 
chosen  as  the  algorithm  of  study.  There  are  several  reasons  for  this  choice.  The  foremost  reason 
is  its  deviation  from  the  standard  Sample  Matrix  Inversion  (SMI)  approach  found  in  Chapter  III. 
By  deviating  from  this  approach,  the  algorithm  brings  into  question  the  applicability  of  work  done 
by  Reed  [15]  on  required  secondary  data  support. 

The  introduction  of  the  CSM  algorithm  begins  with  a  brief  exploration  of  some  of  the  funda¬ 
mental  papers  in  the  field  of  reduced  dimension  STAP.  Prom  here  it  progresses  into  formulating  the 
direct  form  STAP  problem  as  described  in  Chapter  III  into  a  generalized  sidelobe  canceler  model. 
This  model  allows  an  easy  development  of  the  CSM  algorithm. 

4.2  Fundamentals 

In  light  of  these  goals,  let’s  begin  by  mentioning  one  of  the  fundamental  papers  in  this  field. 
This  article  [2]  developed  the  expression  introduced  in  Chapter  III  for  the  optimum  antenna  weights 
in  the  adaptive  processor.  As  shown  by  Brennan  and  Reed,  the  probability  of  false  alarm  Ppx  and 
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probability  of  detection  Pd  for  a  signal  model  given  by  x  =  v  +  n  are 


■Pfa 

Pd 


(4.1) 


where  w  represents  the  complex  valued  weight  vector  for  the  adaptive  processor,  v  the  desired 
signal  vector  or  space-time  steering  vector,  n  the  noise  vector  composed  of  zero  mean  samples  with 
covariance  matrix  R,  and  ^  is  a  normalized  threshold.  The  Q  function  is  defined  as 


Q{a,l3)=  [  ve  (4.2) 

J(3 

where  h{x)  is  the  modified  zero-order  Bessel  function  of  the  first  kind.  These  expressions  were 
developed  by  first  deriving  the  likelihood  ratio  test  for  the  detection  of  a  target  in  the  presence  of 
noise  and  then  applying  the  optimum  Bayes  test.  This  is  the  origin  of  the  parameter  Since  the 
parameter  is  not  of  particular  concern  to  this  discussion,  it  is  only  mentioned  in  passing  here. 

Brennan  and  Reed  then  proceeded  to  maximize  the  above  expression  for  Pd*  They  first 
proved  the  function  given  is  a  monotonically  increasing  function  of  a  =  77^^)  thereby  showing 
the  probability  of  detection  will  be  highest  for  the  filter  w  that  maximizes  a.  The  second  parameter, 
/3,  is  of  no  concern  because  it  is  the  normalized  threshold  from  the  optimum  Bayes  test  and  will 
remain  constant. 

Since  they  had  already  proven  Pd(o;)  is  monotonically  increasing  for  a,  it  follows  the  function 
is  also  increasing  for  Squaring  a  shows  it  is  the  Signal-to-Noise-Ratio  (SNR), 


a 


|w^vp 


w^Rw 


(4.3) 
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The  problem  of  maximizing  the  probability  of  detection  has  now  been  shown  to  be  equivalent  to 
maximizing  the  SNR,  an  important  insight.  Because  the  noise  covariance  matrix  is  Hermitian 
(R  =  R^),  then  it  can  be  diagonalized  through  the  use  of  a  unitary  matrix.  This  is  the  eigenvalue 
decomposition  of  R.  Using  this  decomposition,  Brennan  and  Reed  showed  the  maximum  SNR  is 
obtained  by  letting 

wr^feR^^v*,  (4.4) 


where  fc  is  a  complex  valued  constant.  This  is  the  unconstrained  maximization  of  the  SNR.  The 
maximum  SNR  is  easily  found  through  substituting  (4.4)  into  (4.3), 


max^ 


(4.5) 


Previously  (see  Chapter  III)  the  optimum  adaptive  processor  weights  were  defined  as  R“^v, 
where  the  conjugation  operation  on  the  steering  vector  v  seems  to  be  missing.  The  difference  arises 
because  Brennan  and  Reed  define  the  output  of  the  adaptive  processing  filter  as  2/  =  instead 
of  2/  =  Carrying  this  difference  through  their  equations  results  in  the  same  expression. 

Equation  (4.4)  is  the  basis  for  full  dimension  STAR  Because  of  the  large  dimensionality 
associated  with  the  interference  covariance  matrix  R,  the  computational  complexity  of  inverting 
it  is  still  too  cumbersome  for  real  time  operation.  Furthermore,  the  fiuid  airborne  environment 
necessitates  estimating  R  from  data  obtained  within  a  Coherent  Processing  Interval  (CPI).  The 
estimate  itself  is  difficult  to  obtain  because  the  amount  of  secondary  data  needed  when  the  number 
of  Degrees  Of  Freedom  (DOF)  retains  full  dimensionality. 

Considering  these  practical  road  blocks,  Brennan  and  Reed  discuss  alternative  methods  of 
implementing  the  optimum  adaptive  filter.  In  particular  a  discussion  of  the  method  of  steepest  as- 
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cent  is  proposed.  However,  better  methods  have  arisen  since  the  paper  was  written.  This  discussion 
will  now  focus  on  one  of  those. 


4^3  Translating  Direct  Form  STAP  into  a  Generalized  Sidelobe  Canceler  Form 

Due  to  the  computational  burden  associated  with  retaining  full  dimensionality  and  using 
all  DOF  within  the  adaptive  processor,  an  optimum  method  is  needed  to  generate  a  reduced 
dimension  algorithm.  The  objective  is  to  trade  DOF  for  the  ability  to  estimate  the  covariance 
matrix  and  calculate  its  inverse  in  real  time,  subject  to  both  the  physical  and  computational 
constraints  associated  with  an  airborne  radar  system. 

The  heart  of  a  reduced  dimension  STAP  algorithm  is  throwing  away  DOF.  Determining  which 
DOF  are  the  most  important  for  the  current  interference  plus  noise  covariance  matrix  estimate  is 
the  problem  at  hand.  For  example,  the  Factored  Time-Space  (FTS)  method  simply  throws  away 
all  Doppler  DOF  by  non-adaptively  filtering  in  the  temporal  domain.  Absolutely  no  consideration 
is  given  to  the  impact  this  approach  has  on  the  ability  to  effectively  null  out  the  significant  sources 
of  interference.  Poor  selection  can  result  in  the  adaptive  processor  expending  all  DOF  to  null  out 
a  large  clutter  source,  with  nothing  left  over  to  form  a  mainbeam  in  the  direction  of  the  target. 

Goldstein  and  Reed  have  come  up  with  an  optimum  way  to  select  a  reduced  subset  of  the 
DOF.  This  subset  maximizes  the  ability  of  the  adaptive  processor  to  null  out  large  sources  of 
interference,  i.e.  it  maximizes  the  Signal-to-Interference  plus  Noise  Ratio  (SINR).  First,  begin  by 
defining  the  normalized  space-time  steering  vector  s  as 


(4.6) 


where  v  was  previously  defined  in  Chapter  III  as  the  Kronecker  product  of  the  spatial  and  temporal 
steering  vectors,  a(i?t)  and  respectively. 
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The  purpose  of  any  adaptive  processing  algorithm  for  a  radar  is  to  enhance  target  detection 
through  maximizing  the  probability  of  detection  Po  and  minimizing  the  probability  of  false  alarm 
Pfa*  a  more  complete  discussion  of  these  two  metrics  and  the  associated  hypothesis  testing  is 
found  in  Appendix  C.  Goldstein  breaks  the  detection  problem  into  two  parts. 


•  First,  the  output  SINK  of  the  adaptive  processor  is  maximized  subject  to  the  constraint  of  a 
normalized  response  in  the  direction  of  the  steering  vector. 

•  Second,  the  output  of  the  adaptive  processor  is  then  subject  to  the  hypothesis  testing  mech¬ 
anisms  mentioned  above  and  discussed  in  detail  in  Appendix  C. 

As  developed  by  Reed  [2],  the  unconstrained  optimization  of  the  SINK  was  given  in  Eqn.  (4.4), 
where  the  coefficient  k  was  a  complex  valued  constant.  The  problem  with  this  solution  is  there 
is  incidental  cancellation  of  the  desired  signal  returns.  Applebaum  then  developed  a  method  to 
calculate  the  optimum  adaptive  processor  weights  subject  to  the  constraint  of  a  normalized  response 
in  the  direction  of  the  steering  vector.  These  adaptive  weights  are  given  by  [1] 


R 


(4.7) 


where  the  constraint  is  Wgip^^s  =  1.  Substituting  the  above  expression  into  the  constraint  shows  it 
is  satisfied.  The  output  y  of  the  adaptive  processor  is  then 


y  =  wsinrX  = 


s^R 
s^R-^s ' 


(4.8) 


Next,  direct  form  STAP  is  transformed  into  a  Generalized  Sidelobe  Canceler  (GSC)  form  by 
Goldstein.  Consider  the  use  of  a  unitary  NM  x  NM  transformation  matrix  T,  where  JV  and  M  are 
the  spatial  and  temporal  degrees  of  freedom.  The  resulting  output  SINK  after  this  transformation 
is  identical  to  the  direct  form  processing  architecture,  where  y  =  wU^rX-  Figure  4.1  shows  a  block 
diagram  of  the  direct  form  model. 
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Figure  4.1  The  direct  form  adaptive  processor  block  diagram  showing  the  opti¬ 
mum  constrained  filter  weights. 


The  form  of  the  unitary  transformation  matrix  T  is 

(4.9) 

B 

where  B  represents  a  blocking  matrix.  The  blocking  matrix  projects  the  data  vector  x  the 
null  space  of  the  steering  vector  s,  i.e.  Bs  =  0.  Therefore,  B  blocks  the  desired  signal  and  allows 
all  other  signals  to  pass.  The  size  of  the  blocking  matrix  is  NM  -lx  NM.  Any  full  row  rank 
matrix  satisfying  Bs  =  0  is  a  valid  signal  blocking  matrix  [7].  See  Appendix  B  for  a  discussion  on 
determining  the  blocking  matrix.  Using  these  definitions,  the  transformed  data  vector  is 


(4.10) 


where  b  is  denoted  the  noise-subspace  data  vector  because  it  contains  the  projection  of  the  data  vec¬ 
tor  onto  the  noise  subspace.  A  block  diagram  of  this  new  model  is  shown  in  Fig.  4.2.  Haimovich  [10] 
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X 


Figure  4.2  A  block  diagram  of  the  blocking  matrix  concept  used  in  the  formula¬ 
tion  of  the  GSC  model  of  STAP. 


discusses  determination  of  the  signal  blocking  matrix  for  a  uniformly  spaced  array  and  spatial  steer¬ 
ing  vector. 

The  covariance  matrix  of  the  data  is  found  through  the  following  manipulations, 


R*  =  f  {xx^}  =  £  {TXX^T^}  =  T£  {xx^}  =  TRT".  (4.11) 

Then  substituting  Eqn.  (4.9)  into  Eqn.  (4.11)  results  in 


s^Rs 

s^RB^ 

^hd 

Rx  = 

R 

s 

= 

= 

B 

L  J 

BRs 

BRB^ 

Rb 
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where  this  partitioned  matrix’s  parameters  are  implicitly  defined  when  transitioning  from  the  next 
to  last  step  to  the  result,  i.e. 


cr^  =  s^Rs, 

rbd  =  S{bd*}  =  BRs, 

Rb  =  £{bb^}=BRB^. 

(4.13) 

Now  let  the  transformation  matrix  T  operate  on  the  steering  vector.  Since  the  steering  vector 
was  normalized  previously  a  unit  vector  results, 


s" 

s"s 

Ts  = 

s  = 

B 

Bs 

(4.14) 


Given  these  definitions,  the  optimal  weight  vector  for  the  GSC  adaptive  processing  model 
wgsc  becomes 


wgsc  = 


i^Rr^s 


R^^ei 

efRT^ei 


1 

-Wb 


(4.15) 


where  s  is  the  transformed  steering  vector  and  the  weight  vector  Wb  is  to  be  defined  shortly.  This 
form  of  the  problem  is  best  viewed  in  a  block  diagram  format,  where  the  reasoning  for  the  1  at 
the  beginning  of  the  weight  vector  wgsg  becomes  readily  apparent.  Figure  4.3  shows  the  block 
diagram. 

As  the  figure  visually  illustrates,  the  incoming  signal  is  projected  onto  the  steering  vector. 
This  projection  passes  the  component  of  x  fbat  contains  the  target  and  any  noise  or  interference 


4-8 


that  projects  onto  s.  The  output  of  this  operation  is  the  scalar  containing  the  desired  component, 
d.  Assuming  the  incoming  data  is  the  desired  steering  vector,  i.e.  there  is  a  target  in  the  direction 
the  radar  is  looking,  then  x  becomes  s  and  d  is  equal  to  s^s.  This  is  the  output  of  a  filter  with 
complex  weights  s.  These  weights  give  a  maximum  response  in  the  direction  of  the  steering  vector 
s.  Hence,  the  output  is  an  estimate  of  the  desired  signal. 

The  operation  of  the  GSC  model  can  be  viewed  in  two  steps.  The  upper  branch  estimates  the 
incoming  signal  at  the  direction  of  interest.  Since  there  is  also  noise  and  interference  coming  from 
this  direction  (and  all  directions),  the  target  estimate  is  corrupted  by  the  noise  and  interference. 
The  lower  branch  of  the  GSC  attempts  to  estimate  the  noise  and  interference  from  the  statistics 
of  the  remaining  data  not  oriented  in  the  direction  of  the  steering  vector.  The  incoming  data 
vector  is  projected  onto  the  subspace  of  the  blocking  matrix  B.  This  blocks  any  components  of 
the  incoming  data  vector  matching  the  steering  vector  since  Bs  =  0.  The  weight  vector  Wb  is 
then  applied  to  b,  a  vector  containing  only  the  noise  and  interference.  These  weights  basically 
generate  an  approximation  to  the  interference  contained  in  d  based  on  the  interference  in  b.  This 
interference  estimate  is  then  subtracted  from  the  signal  plus  interference  and  noise  in  the  upper 
branch  resulting  in  an  ideally  pure  target /signal  estimate  corrupted  only  by  noise  applied  to  the 
radar’s  threshold  detection. 

The  optimum  value  of  the  weight  vector  Wb  shown  in  Fig.  4.3  and  Eqn.  (4.15)  is  determined 
by  the  Wiener  solution  [8], 

Wb  =  Rb“^rbd,  (4.16) 

where  the  cross  correlation  vector  and  autocorrelation  matrix  were  defined  previously. 

Figure  4.3  is  the  model  of  fully  adaptive  STAP  using  the  GSC.  Just  as  before,  the  compu¬ 
tational  costs  drive  the  need  for  a  reduced  dimension  application.  Some  of  the  features  of  this 
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model  will  be  used  to  develop  a  reduced  dimension  algorithm  producing  a  maximum  SINK  with 
the  available  DOF. 

Since  the  purpose  of  the  lower  branch  of  the  GSC  is  to  gain  an  estimate  of  the  interference 
contained  in  the  upper  branch,  the  development  now  assumes  x  is  interference  only.  This  allows 
the  analysis  to  parallel  that  of  a  Wiener  filter,  the  upper  branch  being  the  desired  component  d 
(interference  only)  and  the  lower  branch  the  interference  estimate  d.  Using  this  interpretation,  the 
output  of  the  GSC  is  the  error  between  the  interference  estimate  and  the  true  interference  value. 
However,  it  is  important  not  to  lose  sight  of  the  fact  that  the  target  signal  is  generally  included  in 
d  also. 

Care  must  be  taken  to  distinguish  interference  from  noise  in  this  development.  Interference  is 
meant  to  include  all  correlated  undesired  components  while  the  term  noise  refers  to  all  uncorrelated 
undesired  components.  This  distinction  is  made  because  only  the  correlated  components  can  be 
filtered.  This  limitation  is  a  direct  result  of  the  filter  being  based  on  the  statistics,  or  correlation,  of 
the  secondary  data.  For  the  interference  estimate  to  accurately  refiect  the  true  interference  in  the 
desired  signal,  the  interference  samples  in  the  secondary  data  must  be  correlated  to  the  interference 
found  in  the  desired  signal. 

Now  that  a  suitable  model  is  available,  let’s  analyze  it’s  use  in  a  reduced  dimension  application. 
From  this  analysis,  we  can  then  choose  a  suitable  metric  for  selecting  which  DOF  maximize  the 
output  SINK  with  reduced  dimensionality.  First,  express  the  covariance  matrix  Rb  for  the  data 
vector  b  in  terms  of  its  eigenvalue/eigenvector  decomposition.  This  decomposition  is  given  by 

Rb  =  UAU^,  (4.17) 

where  U  is  the  unitary  matrix  containing  the  eigenvectors  and  A  is  the  diagonal  matrix  containing 
the  eigenvalues.  In  the  derivation  of  the  following  equations,  the  unitary  property  allows  the  use  of 
U^U  =  UU"  =  I.  The  size  of  U  and  A  is  {NM  —  1)  x  {NM  —  1).  Now,  if  we  project  the  vector 
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b  to  a  different  coordinate  space  using  the  unitary  matrix  of  U^, 


p  =  U^b,  (4.18) 

then  we  can  start  viewing  the  matrix  U  in  terms  of  reducing  the  dimension  of  the  problem.  Let’s 
first  examine  the  statistics  of  p,  the  vector  Goldstein  [6]  refers  to  as  a  principal  coordinate  process. 
The  covariance  matrix  Rp  of  p  is  given  by 


Rp  =  5  {pp^}  =  U^RbU  =  A. 


(4.19) 


The  cross  correlation  between  p  and  d  is 


rpd  =  f{pd*}  =  U"rbd.  (4.20) 

And  finally,  the  Wiener-Hopf  solution  for  the  optimum  weight  vector  Wp  in  the  full  dimension  GSC 
with  the  eigenvector  matrix  projection  is 


Wp  =  Rp  ^rprf  =  U^Wb. 


(4.21) 


A  block  diagram  of  the  full  dimension  GSC  model  using  the  full  eigenvector  matrix  as  a 
projection  is  shown  in  Fig.  4.4.  The  output  of  the  GSC  model  can  be  inferred  just  from  the 
diagram  alone  or  by  going  back  through  the  equations.  This  output  is 

2/  =  (s^-w^B)x.  (4.22) 

At  this  point  in  Goldstein’s  derivation,  the  GSC  model  is  still  using  the  full  dimensionality 
of  the  problem.  DOF  reduction  is  performed  through  the  use  of  the  unitary  transformation  matrix 
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X 


Transformation  Matrix,  T 

Figure  4.4  The  full  dimension  GSC  block  diagram  with  eigenvector  projection. 


of  eigenvalues  U.  In  the  case  of  reduced  dimensionality,  only  a  subset  of  the  columns  of  U  are 
used.  This  subset  of  the  total  number  of  eigenvectors  gives  the  new  transformation  U  a  size 
of  (NM  -  1)  X  0,  where  U  is  composed  of  O  columns  of  U.  Multiplying  b  by  M  reduces  the 
dimension  to  O.  A  block  diagram  of  the  reduced  dimension  GSC  model  is  shown  in  Fig.  4.5.  As 
inferred  by  the  figure,  the  output  of  the  dimension  reduction  z  is  defined  as 


(4.23) 


The  next  problem  is  determining  which  eigenvectors  represent  the  best  choice  for  the  con¬ 
struction  of  W.  Goldstein  begins  this  process  by  looking  at  the  reduced  dimension  covariance  matrix 
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Transformation  Matrix,  T  Dimension  Reduction 

Figure  4.5  The  reduced  dimension  GSC  block  diagram  using  O  eigenvectors  as  a  projec¬ 
tion. 
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R,  =  W^RbW  =  Ao, 


(4.24) 


Rz) 


where  Ao  is  the  diagonal  matrix  composed  of  the  eigenvalues  corresponding  to  the  eigenvectors 
used  to  construct  ll.  The  cross  correlation  vector  between  the  reduced  dimension  noise  subspace 
vector  z  and  d  is 


Tad  =  S  {zd*}  = 


(4.25) 


Using  these  statistics,  the  Wiener  solution  w*  for  the  O  dimensional  weight  vector  is 


Wgf  =  Rz  ^rzd  =  Ao^W^rbd.  (4.26) 

All  of  the  structures  necessary  to  build  the  adaptive  filter  are  in  place.  The  only  issue  in  need 
of  consideration  is  the  choice  of  eigenvectors  used  to  compose  the  dimension  reducing  matrix  U. 
Goldstein  states  [8]  the  most  popular  technique  for  this  subspace  selection  is  based  on  the  principal 
components  method. 

The  principal  components  method  chooses  the  set  of  eigenvectors  resulting  in  the  closest 
approximation  to  the  the  full  covariance  matrix  Rb-  This  is  the  set  of  eigenvectors  corresponding 
to  the  O  largest  eigenvalues  of  the  matrix.  The  problem  with  this  method  is  that  it  results  in  the 
optimum  representation  of  the  covariance  matrix.  This  is  not  the  point  of  the  problem. 

The  eigenvectors  should  be  chosen  to  maximize  the  SINK  of  the  adaptive  filter.  As  shown  by 
Goldstein,  this  selection  is  not  necessarily  the  same  as  choosing  the  eigenvectors  corresponding  to 
the  largest  eigenvalues.  This  selection  is  more  than  just  a  function  of  the  covariance  matrix,  it  is 
also  a  function  of  the  cross  correlation  between  the  beamformed  signal  d  and  the  noise  process  b. 
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Goldstein  provides  the  derivation  of  the  optimum  selection  of  the  eigenvectors  in  another 
paper  [6].  The  eigenvectors  are  selected  to  satisfy  a  Minimum  Mean  Square  Error  (MMSE)  criterion. 
The  output  of  the  full  dimension  GSC  model  was  previously  shown  in  Eqn.  (4.22)  to  be 

l,  =  (s^-w^B)x  (4.27) 

and  the  Wiener-Hopf  solution  was  shown  in  Eqn.  (4.16)  as 


Wb  =  Rb  ^rbd. 


(4.28) 


Using  these  two  equations,  the  MMSE  P  is 

P  =  ^  {IJ/P}  =  Od-  rbrfR-b'^rbd,  (4.29) 

where  the  cross  correlation  vector  and  covariance  matrix  of  the  transformed  data  vector  were 
defined  previously.  The  MMSE  is  conserved  by  a  unitary  transformation  [6],  therefore  the  MMSE 
associated  with  the  full  dimension  GSC  model  shown  in  Fig.  4.4  remains  the  same. 

Goldstein  uses  the  term  MMSE  because  we  are  considering  the  case  of  noise  interference  only. 
Using  this  assumption,  the  output  of  the  GSC  model  should  be  zero  given  a  perfect  estimate  of 
the  interference.  Hence,  at  this  point  in  the  derivation  any  output  of  the  GSC  is  an  error  in  the 
interference  estimate. 

Denote  the  output  of  the  reduced  dimension  GSC  model  as  2/rj  following  Fig.  4.5.  The  output 
in  terms  of  the  transformed  data  vector  b  is 


=  d-y^—  d  - 


(4.30) 
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Adding  and  subtracting  the  output  of  the  full  dimension  GSC  does  not  modify  the  equation, 


yr  =  y-y  +  d-yz. 


(4.31) 


Referencing  Fig.  4.4  shows  we  can  write  2/  as  d  -  2/p.  Therefore,  substituting  for  one  of  the  terms 
gives 


yr  =  y-d  +  yp  +  d-y„-y-hyp-yz, 


(4.32) 


or  in  terms  of  the  matrices  already  defined, 


r  n 

y 

yr  = 

b 

(4.33) 


Using  this  expression,  the  {NM  —  1)  x  1  vector  e  representing  the  error  between  the  reduced 
dimension  weight  vector  and  the  full  dimension  weight  vector  is  defined  by  Goldstein  as 


e  =  Uwp  -  Wwa. 


(4.34) 


With  these  forms,  the  mean-square  value  of  the  reduced  dimension  GSC  model  output  is 


^  {Ij/rD  =£’{|y  +  e^bp}. 


(4.35) 
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Expanding  the  quadratic  uncovers  the  error  found  in  the  full  dimension  GSC  model  plus  terms  due 
to  the  reduced  dimensionality, 

+  5  {(2/  +  e^b)(j/  +  e^b)^} 

=  S^\y\^  +  yh^e  +  y*e^h  +  e^bb^ej 

=  P  +  e^Rb^e  +  S  {j/b^e}  +  €  {y*e^b}  .  (4.36) 

Ridding  ourselves  of  the  last  two  terms  is  interesting.  Using  Fig.  4.4  as  a  reference,  examine  the 
purpose  of  the  model.  The  only  input  under  consideration  at  this  point  is  that  of  the  interfernce. 
Therefore,  the  lower  branch  is  an  estimate  of  the  interference  in  the  upper  branch  and  the  output 
y  is  simply  the  difference  between  the  estimate  and  the  true  value, 

y  =  d-l  (4.37) 

Again  using  the  figure,  the  estimator  is  based  on  the  noise/interference  projection  vector  b  such 
that  d  =  w^TJ-^b.  The  full  dimension  Wiener  filter  was  already  determined  in  Eqn.  (4.21)  as 

Wp  =  Rp“^rprf  =  U^Wb.  (4.38) 


Therefore,  the  noise  estimator  becomes 

d  =  w^UU^b  =  w^b. 


(4.39) 
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where  the  property  UU^  =  I  is  invoked  for  the  unitary  matrices.  Using  the  Wiener  solution  for 
the  above  filter  found  in  Eqn.  (4.16),  the  estimator  can  be  rewritten  as 

d  =  r^^Rb-'b,  (4.40) 

where  the  Hermitian  property  of  the  covariance  matrix  has  been  invoked,  i.e.  Rb^  =  Rb- 

Now  let’s  examine  the  third  term  of  Eqn.  (4.36)  using  the  expression  for  the  estimate  of  d, 

£  {j/b^e}  =  5  I  b^l  e 

=  5{(d-r^,Rb-^b)b^}e 
=  £:{db^}e-£:{rg,Rb“^bb^}e 

=  rdbC  -  rbdRb~^Rbe 
=  rb<ie-r£,e 

=  0,  (4.41) 


where  the  final  result  uses  the  fact  r^b  =  The  fourth  term  in  Eqn.  (4.36)  can  be  shown  to  be 
equal  to  the  Hermitian  of  the  term  just  shown  to  equal  zero, 

S  =  S  {yh^e}^  ,  (4.42) 

hence  it  is  also  zero. 

It  is  interesting  to  note  the  above  mathematics  were  not  necessary  to  show  the  two  terms  are 
zero.  Since  the  estimate  of  d  is  the  result  of  a  Wiener  filtering  operation,  then  by  definition  the 
error  is  minimized.  This  minimization  is  another  way  of  stating  the  error,  y  =  d  —  d,  is  orthogonal 
to  the  data  vector  b  used  to  generate  the  statistics  for  the  Wiener  filter.  Since  y  at  this  point  is 
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an  expression  of  the  error  between  d  and  d,  then  y  is  orthogonal  to  b  and  their  inner  product  is 
obviously  zero. 

Either  through  the  mathematical  manipulation  or  the  Wiener  filtering  rationalization,  the 
expression  in  Eqn.  (4.36)  reduces  to  the  one  given  by  Goldstein  [6] 

^  {l2/rp}  =s{\y  +  e^b|^}  =  P  +  e^RbC,  (4.43) 

where  P  is  the  mean  square  output  of  the  full  dimension  GSC  model.  This  implies  the  error  due 
to  the  dimension  reduction  is  contained  in  the  second  term,  e^RbC.  Minimization  of  the  error  is 
done  via  minimization  of  this  term. 

At  this  point  in  the  derivation,  we  deviate  from  Goldstein’s  development  to  one  of  our  own. 
Let  us  begin  the  minimization  by  first  expressing  the  error  term  explicitly, 

e^Rbe  =  (w^ U"  -  wf  W")  Rb  (Uwp  -  Ww, ) .  (4.44) 

The  covariance  matrix  Rb  was  the  subject  of  the  previous  eigenvalue/eigenvector  decomposition. 
Using  this  decomposition,  we  can  substitute  into  Eqn.  (4.44), 

e^Rbe  =  U"  -  UAU"  (Uwp  -  Ww.) .  (4.45) 

Distributing  the  unitary  eigenvector  matrices  results  in 

e^Rbe  =  (w" U^U  -  wf  M^u)  A  (U^Uwp  -  U"Ww*)  (4.46) 
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and  allows  us  to  invoke  the  unitary  property  of  the  eigenvector  matrix  U  to  further  simplify  the 
error  expression  to 

e^Rbe  =  (w"  -  A  (wp  -  U^Mw*) .  (4.47) 

If  we  pull  the  hermitian  operation  out  so  it  operates  on  the  entire  left  quantity  and  factor  the 
diagonal  eigenvalue  matrix  we  see  the  equation  is  basically  the  magnitude  squared  of  a  simpler 
term, 


e^Rbe  =  (wp  -  U^Uwz)^  A2  (wp  —  U^Wwz) 

=  A2  (wp  -  U^Uwg)  I  , 


(4.48) 


where  the  |  |  operator  denotes  the  magnitude  of  the  vector.  Because  the  covariance  matrix  R  is 

hermitian,  then  Rb  is  also  hermitian.  For  any  hermitian  matrix,  the  eigenvalues  are  real  valued  [3, 

page  119].  This  allows  us  to  take  the  hermitian  of  the  diagonal  eigenvalue  matrix  A^  without  any 
r  1 

effect,  i.e.  \  =  A2 .  Using  Eqn.  4.48,  minimization  of  the  error  term  becomes 


min  (e^Rbc)  =  min  |  A 2  (wp  -  U^Ww*)  |  =  min  |  A 2  (wp  -  U^Wwz) || .  (4.49) 


Now  let’s  begin  the  analysis  of  this  equivalent  minimization  term  in  the  hopes  of  further 
simplifying  it  and  developing  an  optimum  method  to  choose  the  eigenvectors  comprising  the  matrix 
U.  First,  distribute  the  eigenvalue  matrix  across  the  quantity, 


min  (e^Rbc)  =  min  |  A  Wp  -  AW^Mwzj  | , 


(4.50) 
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Prom  Eqn.  (4.21)  we  can  substitute  for  Wp, 


min  (e^Rbe)  =  min  | 


A2U^Wb-A2U^Mw* 


(4.51) 


Using  Eqn.  (4.16)  to  substitute  for  the  NM  —  1  dimensioned  weight  vector  Wb  and  Eqn.  (4.26)  for 
the  0  dimensional  weight  vector  results  in 


min  (e^Rbe)  =  min  1 1 A^ U^Rb  ^rbd  -  A^ U^WAo^W^rbrf|  | . 


(4.52) 


Since  we  previously  defined  the  eigenvector /eigenvalue  decomposition  of  Rb  as  UAU^,  then  the 
inverse  of  Rb  is 


Rb“^  =  UA-^U^, 


(4.53) 


allowing  us  to  substitute  into  Eqn.  (4.52)  and  again  invoke  the  property  of  the  unitary  eigenvector 
matrix, 


min  (e^Rbe) 


min  ||A=U^UA-^U"rbd  -  A^U^WAo^W^rw 
min  ( |A-  =  U"rb<i  -  A^U^WA^^M^rbrifl . 


(4.54) 


The  next  step  in  the  development  involves  the  ordering  of  the  eigenvectors  and  eigenvalues 
in  the  decomposition  of  Rb.  This  decomposition  can  be  expressed  in  the  form  of  a  summation, 

JVM-l 

Rb  =  UAU^=  Xmvf.  (4.55) 

2=0 


Examining  this  form  shows  the  ordering  of  the  columns  of  U  and  A  is  totally  insignificant  as  long  as 
the  eigenvalue/eigenvector  relationship  is  maintained.  At  this  point  in  the  derivation,  let  us  assume 
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the  decomposition  is  ordered  such  that  the  first  O  columns  of  U  and  A  contain  the  eigenvectors 
and  eigenvalues  corresponding  to  li.  Then  these  two  matrices  can  be  written  in  the  form, 


U 


u  c 


(4.56) 


A 


Ao  0 
0  Anm-1-o 


(4.57) 


where  C  is  the  matrix  of  eigenvectors  not  contained  in  W,  Anm-i-o  is  the  diagonal  matrix  of 
eigenvalues  corresponding  to  C,  and  0  represents  a  matrix  of  zeros. 

Substituting  this  form  of  U  into  Eqn.  (4.54)  allows  us  to  write 


min  (e^Rbe) 


=  mm 


=  mm  < 


A-^UVd-A* 


A-iu^Fbd-Ai 


- 

2n 

► 

> 

U^U 

c^u 


Ao^W^rbd 


r 

■ 

2  > 

=  min  < 

A-»U"rbd-A^ 

lo 

Ao^W^rbd 

1 

0 

j 

(4.58) 


where  the  final  step  takes  advantage  of  the  orthogonality  of  the  eigenvectors.  Furthermore,  multi¬ 
plying  the  matrix  [I©  ;  0]  by  simply  picks  off  the  first  O  eigenvalues, 


min  (e^Rbe)  =  min  < 


A-*U^rb<i  - 


Ao'MVd 


(4.59) 
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allowing  us  to  further  simplify  the  minimization  by  performing  one  more  matrix  multiplication, 


min  (e^Rbe)  =  min  < 


f 

- 

A-^U^Pbd  - 

Ao" 

U^Tbd 

1 

k 

0 

J 

(4.60) 


To  find  the  optimum  method  for  selecting  the  eigenvectors  to  compose  W,  we  need  to  expand 
the  vector  magnitude  operation  in  Eqn,  (4.60)  to  give  a  minimum  error,  min  (e^Rbe) ,  equal  to 


mm  < 


< 

/ 

\ 

H 

A-  =  U^rbd  - 

O  1 

U^Tbd 

A-»U"rbd  - 

Ao’ 

U^Thd 

> 

k 

1 

0 

0 

) 

J 

(4.61) 


Distributing  the  hermitian  operation  throughout  the  first  quantity  and  using  the  fact  that  the 
eigenvalues  are  real^  results  in  the  error  expression  taking  the  form 


f 

y 

min  < 

Ao"  ol) 

A-^U^rbd  - 

iC4 

O 

<! 

U^T^d 

\ 

k 

L  ^  J/  1 

1 

[»J 

) 

J 

(4.62) 


Further  insight  into  the  algorithm’s  error  requires  multiplying  the  two  quantities  as  shown  by 


mm 


in  (e^Rbc)  =  min  < 


r£jUA-5A-5U®rbd  +  r£,W 


Aq*  0 


Ao" 

0 


U^Thd 


^bdU 

O  1 

o 

A-  =  UVd-rg,UA-i 

O  1 

L  J 

[ » J 

J 

^The  eigenvalues  are  real  because  the  eigenvalue  decomposition  was  performed  on  the  hermitian  matrix  Rb- 
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Reducing  this  equation  is  a  matter  of  recognizing  we  can  substitute  the  vector  magnitude,  |  |^,  of 
each  quantity.  First,  rearrange  Eqn.  (4.63)  so  that  min  (e^Rbe)  is  equal  to 


_1 

H 

_1 

\ 

min  \ 

(A-^U^rbd)"  (A-^U'^rbrf)  -|- 

U^T^bd 

Ao^ 

U^Vbd 

[ 

0 

) 

0 

/ 

(4.64) 


The  first  two  terms  are  recognized  as  the  magnitude  of  a  vector.  The  second  two  terms  require 
some  manipulation  to  see  they  are  also  vector  magnitudes.  The  following  shows  the  first  step  in 
the  manipulation  of  these  terms, 


min  (e^Rbe)  =  min||A  ^u^rbd 


-^bd^  aJ  0 


+ 


Ao  =  M^rbd| 

U^rbd  -  r£,U 


Ao" 


-i 


Ao'M^rbd 


}  .  (4.65) 


Distributing  the  eigenvector  matrices  U  and  14  results  in 


min(e^Rbe)  =  min||A  ^U®rb<i|  +  |Ao  =  W-®^rb<i 

0 


-^bd  MAo“  0 


rbd  -  v^JUAq  ^  Aq  ’  U^v\>d 


(4.66) 
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where  we  have  taken  advantage  of  the  fact  that  U  was  previously  defined  to  have  the  form  of 
U  =  [  W  C  ].  This  form  of  the  equation  shows  the  last  two  terms  are  also  vector  magnitudes, 


min  (e^Rbc)  =  min  (  A  ^U^rbd 

2 

+ 

Aq  *  U^Tbd 

[ao"H 

rbd 

0 

2 

A-Q^U^Tbd 

(4.67) 


Since  the  second  and  fourth  terms  cancel,  the  next  operation  necessary  is  to  simply  distribute  the 
cross  correlation  vector  rbd  on  the  third  term  as  shown  by 


min  (e^Rbc)  =  min  < 


A-iu^rbrf 


(4.68) 


The  second  term  simplifies  even  further  because  it  is  the  magnitude  of  a  vector.  The  0  matrix  does 
not  affect  the  magnitude  of  the  vector  and  can  be  dropped  resulting  in 


min  (e^Rbc)  =  min  |  j  A  2  U^rbd 


(4.69) 


The  final  step  merely  substitutes  the  summation  form  of  the  vector  products  shown  as  the  operands 
of  the  I  operation, 


min  (e^Rbc)  =  min  < 


NM-l  jT 


i-1 


(4.70) 


Prom  this  equation,  the  O  eigenvectors  comprising  U  should  be  selected  according  to  the  metric 


\/Aj 


(4.71) 
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Choosing  the  O  eigenvectors  corresponding  to  the  0  maximum  values  of  this  metric  will  produce 
the  minimum  error  in  the  reduced  dimension  algorithm,  thereby  producing  the  highest  SINR. 

This  is  the  same  result  as  shown  by  Goldstein’s  work  [6],  which  he  designates  the  cross  spectral 
metric.  Due  to  the  use  of  this  metric,  the  algorithm  is  referred  to  as  the  Cross  Spectral  Metric 
(CSM)  algorithm  for  reduced  dimension  STAR  Obviously,  this  greatly  differs  from  the  principal 
components  method.  If  the  eigenvectors  of  Rb  are  rank  ordered  vector  according  to  Eqn.  (4.71), 
then  the  selection  of  U  is  simply  a  matter  of  selecting  the  first  0  eigenvectors. 

The  correlation  between  the  desired  signal  and  the  noise  projection  vector  T^d  is  interpreted 
as  the  cross-spectral  energy.  Therefore,  the  metric  basically  measures  the  cross-spectral  energy 
projected  along  the  eigenvector.  This  energy  relationship  explains  the  squaring  of  the  metric 
as  shown  in  Eqn.  (4.71).  This  energy  comparison  is  an  important  consideration  when  trying  to 
understand  the  advantages  of  the  algorithm.  The  selection  process  solves  the  problem  of  inadver¬ 
tently  choosing  an  eigenvector  along  which  the  signal  has  no  energy,  i.e.  the  signal  is  naturally 
experiencing  a  null  in  this  direction.  Because  there  is  a  natural  null,  the  noise  has  no  effect  on 
the  signal.  Expending  valuable  DOF  creating  a  null  in  the  antenna  pattern  to  filter  this  noise  is  a 
waste.  Creating  antenna  pattern  nulls  in  directions  of  natural  nulls  is  the  problem  inherent  in  the 
principal  components  method  of  eigenvector  selection. 

As  validation  of  the  mathematics.  Fig.  4.6  compares  the  SINR  for  both  methods  using  the 
interference  plus  noise  covariance  matrix  R  from  the  previous  examples.  The  expression  for  the 
SINR  is  developed  in  Chapter  VI,  see  Eqn.  (6.20).  The  figure  is  included  here  to  illustrate  the 
improvement  in  algorithm  performance  when  the  CSM  is  used. 

Figure  4.6  shows  the  CSM  method  of  choosing  U  clearly  outperforms  the  method  of  principle 
components.  The  algorithm  rapidly  converges  to  the  optimum  fully  adaptive  SINR  and  outperforms 
the  principal  components  method  at  all  dimensions  up  to  the  point  where  both  reach  the  optimum 
SINR. 
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SINR  Comparison 


Figure  4.6  A  comparison  of  the  reduced  dimensionality  SINR  for  the  principal 
components  and  cross  spectral  metric  methods.  The  interference  plus  noise  covari¬ 
ance  matrix  R  is  the  same  as  the  one  used  to  generate  the  SINR  plots  for  Section  3.5. 
The  cross  spectral  method  results  in  a  more  rapid  convergence  to  the  optimum  (fully 
adaptive)  SINR  than  the  principal  components  method. 
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SINRcsm  part  (dB) 


In  summary,  Goldstein  has  shown  the  principal  components  method  does  not  produce  the 
optimum  SINK  for  a  reduced  dimension  adaptive  processor.  Instead,  the  GSM  model  as  derived 
by  Goldstein  results  in  the  maximum  SINK. 

Summary 

In  order  to  gain  insight  into  the  problem  of  choosing  the  best  DOF  for  reduced  dimensionality 
STAP,  the  GSM  model  by  Goldstein  was  examined  in  detail.  This  led  to  the  cross-spectral  metric, 
a  method  that  results  in  the  maximum  SINK  for  a  reduced  dimension  algorithm. 

The  basis  for  the  cross  spectral  metric  is  the  use  of  a  subset  of  the  eigenvectors  of  the 
interference  plus  noise  covariance  matrix  Rb  as  the  reduced  dimension  projection  matrix.  The 
cross  spectral  metric  then  produces  the  optimum  SINK  for  this  projection  matrix.  If  a  different 
projection  matrix  is  used,  it  is  reasonable  the  SINK  could  be  further  improved. 

The  GSM  algorithm  is  computationally  expensive.  The  NM  —  1  dimension  covariance  matrix 
Rb  must  first  be  determined.  Further  compounding  the  burden  is  the  eigenvalue/eigenvector 
decomposition  of  this  matrix,  required  to  build  the  reduced  dimension  projection  matrix  U,  Savings 
in  required  secondary  data  support  are  not  obvious  because  they  come  about  through  the  Maximum 
Likelihood  Estimator  (MLE)  used  to  estimate  the  covariance  matrix  R. 

If  the  secondary  data  support  required  for  the  GSM  algorithm  follows  the  work  done  by 
Reed  [15],  near  optimum  SINK  performance  should  be  obtained  by  only  using  20  secondary  data 
vectors  for  the  estimate  of  R.  Without  the  reduced  dimensionality,  NM  secondary  data  vectors 
are  required  to  reach  optimum  performance.  Ghapter  V  delves  deeper  into  the  MLE  for  R  and  also 
introduces  the  problems  associated  with  non-homogeneities  in  the  secondary  data  used  to  estimate 
R.  Furthermore,  a  method  of  excising  potential  non-homogeneities  is  discussed  and  analyzed  in 
detail  for  the  FTS  method.  Ghapter  VI  experimentally  examines  the  GSM  algorithm  in  terms  of 
required  secondary  data  support  and  the  impact  of  non-homogeneities  in  the  data. 
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V.  Covariance  Matrix  Estimation  and  Non-homogeneity  Detection 

5  A  Introduction 

One  of  the  primary  difficulties  encountered  in  the  implementation  of  a  practical  adaptive 
processor  is  estimation  of  the  noise  covariance  matrix.  All  of  the  algorithms  discussed  rely  on  the 
interference  plus  noise  covariance  matrix  R  for  the  adaptive  weight  calculation.  Because  of  the 
rapidly  changing  environment  encountered  in  an  airborne  platform,  the  noise  statistics  must  be 
estimated  from  the  data  available.  Furthermore,  this  estimate  must  be  repeated  for  each  range  cell 
of  interest.  This  repetition  is  necessary  because  the  underlying  assumption  when  examining  each 
range  cell  is  that  there  is  a  potential  target  within  the  range  cell.  This  potential  target  represents 
a  non-homogeneity  in  the  data  and  should  not  be  part  of  the  estimate,  even  though  a  target  does 
not  necessarily  exist  at  the  range  gate  of  interest. 

The  data  can  be  broken  up  into  two  main  classes,  the  primary  data  consisting  of  the  range 
gate/cell  of  interest  and  the  secondary  data  consisting  of  everything  else.  The  secondary  data 
is  used  for  estimating  the  noise  covariance  matrix.  Anything  other  than  noise  represents  a  non¬ 
homogeneity  in  the  data.  These  non-homogeneities  can  range  in  severity  from  a  surface  clutter 
transition  (e.g.  land  clutter  to  sea  clutter)  to  a  strong  target  in  a  range  cell  other  than  the  one  of 
interest. 

In  the  adaptive  filter  developed  to  this  point  using  the  CSM  algorithm  or  the  development 
shown  in  Chapter  III,  the  primary  goal  is  to  filter  out  the  sources  of  interference.  Any  sources 
contained  in  the  estimated  covariance  matrix  will  result  in  the  adaptive  filter  expending  Degrees 
Of  Freedom  (DOF)  to  build  an  antenna  pattern  with  nulls  on  these  sources.  Some  examples  of 
interference  sources  desired  in  the  covariance  matrix  include  strong  sources  of  clutter  and  barrage 
noise  jamming.  Conversely,  the  estimate  should  not  include  other  targets.  These  targets  waste 
valuable  DOF  building  nulls  with  no  real  advantage  in  the  detection  of  targets  in  the  direction  of 
interest. 
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5.2  Covariance  Matrix  Estimation 


The  interference  plus  noise  covariance  matrix  R  is  the  foundation  of  the  Wiener  filter  solution 
for  the  optimum  weight  vector  w.  In  Chapter  III,  this  matrix  was  assumed  known.  In  the  exam¬ 
ples  shown  in  that  chapter,  it  was  estimated  using  the  models  developed  in  Sections  3.3.1,  3.3.2, 
and  3.3.3.  However,  the  radar  is  now  placed  in  the  real  world  environment  where  all  the  parameters 
are  constantly  changing.  The  interference  plus  noise  covariance  matrix  must  be  estimated  from  the 
data  available. 

The  Maximum  Likelihood  Estimate  (MLE)  R  for  the  interference  plus  noise  covariance  matrix 
is 

/+J/2 

R«  =  4  E  XiXf  (5.1) 

i=l-J/2 

where  the  I  subscript  explicitly  shows  the  estimate  changes  for  each  range  gate,  J  is  the  number  of 
secondary  data  vectors  to  be  used  in  the  estimate,  and  Xi  is  the  i^^  range  gate  space-time  snapshot 
or  data  vector.  As  done  previously  in  Chapter  III,  the  dependence  on  range  will  be  implicit  from 
this  point  forward  and  the  I  subscript  dropped.  Obviously,  the  best  estimate  is  obtained  using  the 
maximum  amount  of  secondary  data  vectors.  This  maximum  is  equal  to  the  number  of  range  gates, 
L. 

The  interference  covariance  matrix  R  in  the  Wiener  solution  or  weight  vector  for  all  of  the 
STAP  algorithms  is  replaced  by  the  estimate,  R.  This  replacement  for  fully  adaptive  STAP  is 

W  =  (5.2) 

where  gt  is  the  tapered  space-time  steering  vector  to  the  target.  If  no  taper  is  desired,  vt  is 
substituted. 
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An  expression  for  the  Signal-to-Interference-plus-Noise  Ratio  (SINK)  is  obtained  by  substi¬ 


tuting  Eqn.  (5.2)  into  Eqn.  (3.91)  to  obtain 


smR  = 


gfR-V* 


gfR-iRR-igt 


(5.3) 


The  problem  with  the  use  of  this  expression  is  the  true  interference  covariance  matrix  R  is  generally 
not  known.  This  is  the  reason  for  using  the  MLE  in  the  first  place.  However,  an  estimate  for  the 
true  interference  covariance  matrix  R  can  be  obtained  by  using  all  of  the  data  vectors  available  to 
the  MLE,  whereas  R  would  be  obtained  with  some  subset  of  these  secondary  data  vectors. 

The  obvious  question  arising  from  the  use  of  this  MLE  concerns  the  performance  degradation. 
This  performance  degradation  is  directly  related  to  the  number  of  secondary  data  vectors,  J,  used 
in  the  MLE.  Reed  found  this  loss  in  SNR  [15]  as 


isNR  =  -lOlogio 


J  +  2  “  iVboF 
J  +  l 


(5.4) 


where  iVooF  is  the  number  of  Degrees  Of  Freedom  (DOF).  Using  this  equation,  performance  within 
3  dB  of  the  optimal  SNR  given  a  known  covariance  matrix  can  be  obtained  by  choosing  J  equal 
to  [15] 


J  =  2iVDOF  —  3  2iVDOF-  (5«5) 

This  is  equivalent  to  a  total  secondary  data  set  size  of  approximately  2iVDOF*  However,  this  result 
was  developed  by  Reed  for  Sample  Matrix  Inversion  (SMI)  algorithms,  i.e.  the  fully  adaptive 
algorithm  discussed  in  Chapter  III.  The  nature  of  the  FTS  reduced  dimension  algorithm  also 
allows  the  use  of  Eqn.  5.5.  The  particular  point  allowing  the  use  of  this  sample  size  approximation 
is  the  FTS  method  computes  the  covariance  matrix  after  the  Doppler  filtering  operation.  This 
filtering  operation  simply  results  in  a  new  data  set  for  the  adaptive  algorithm.  From  that  point 
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onward,  the  computation  of  the  weight  vector  is  similar  to  that  of  a  fully  adaptive  algorithm.  The 
interference  covariance  matrix  is  estimated  from  the  Doppler  filtered  data  and  then  inverted  in  the 
process  of  forming  the  weight  vector.  The  technique  still  falls  under  the  SMI  umbrella. 

The  CSM  algorithm  developed  by  Goldstein  differs  from  the  SMI  techniques  introduced  by 
Ward  [21]  and  Reed  [2,15].  As  discussed  in  Section  4.3,  the  reduced  dimension  CSM  algorithm 
requires  the  full  dimension  covariance  matrix.  The  eigenvalue  decomposition  of  this  matrix  is 
performed  and  then  used  to  generate  a  reduced  dimension  covariance  matrix.  This  matrix  is  then 
inverted  in  the  computation  of  the  optimum  weight  vector.  The  concerns  with  applying  Eqn.  5.5  to 
the  CSM  reduced  dimension  algorithm  are  analyzed  in  detail  in  Chapter  VI.  This  chapter  focuses 
on  understanding  the  MLE  and  its  limitations. 

A  more  fundamental  problem  with  the  MLE  estimator  is  its  assumption  of  independent, 
identically  distributed  (iid)  data.  This  assumption  is  not  generally  valid.  This  is  the  reason  for  the 
non-homogeneity  detector  proposed  by  Melvin  [14]. 

5.3  Non-Homogeneity  Detection 

As  just  explained,  the  MLE  is  based  on  the  assumption  the  secondary  data  is  iid.  Generally, 
this  assumption  is  not  valid  for  real  world  data.  Therefore,  Melvin  [14]  introduced  a  method  to 
select  the  vectors  out  of  the  entire  set  of  secondary  data  that  conform  best  to  the  iid  assumption. 

Melvin  begins  by  defining  homogeneous  data.  Given  the  covariance  estimate  from  two  differ¬ 
ent  range  cells  i  and  j, 


=  £{XiXi}  and  Rj  =  £{XjXf}, 


(5.6) 
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then  a  homogeneous  set  of  data  satisfies  [14] 


RiR^.  ^  «  Inm^  (5-7) 

where  NM  x  NM  is  the  dimension  of  the  interference  covariance  matrix.  The  number  of  pulses 
in  a  Coherent  Processing  Interval  (CPI)  and  number  of  antenna  channels  are  denoted  M  and  AT, 
respectively. 

Considering  the  other  side  of  the  problem,  the  secondary  data  vector  is  defined  in  terms  of 
the  components  it  contains, 

Xi  Xci  d”  Xji  'X.ni  Xoi*  (^•^) 

The  components  representing  the  homogeneous  data  are  the  homogeneous  clutter  component  Xch 
white  noise  jamming  Xji^  uncorrelated  noise  Xni-  The  non-homogeneous  component  is  repre¬ 
sented  by  Xoi  consists  of  colored  noise  jamming  and  the  non-homogeneous  clutter  component. 

Melvin  proposes  the  Generalized  Inner  Product  (GIP) 

9i  =  (5.9) 

as  a  method  to  determine  a  secondary  data  vector’s  degree  of  homogeneity.  The  matrix  R^  repre¬ 
sents  an  interference  plus  noise  covariance  matrix  estimate  based  on  the  entire  secondary  data  set. 
The  GIP  is  then  calculated  for  every  vector  contained  within  the  secondary  data  set,  producing  the 
scalar  gi.  The  homogeneous  data  to  be  used  in  calculating  the  new  interference  covariance  matrix 
estimate  is  then  determined  by  choosing  those  vectors  corresponding  to 


g  -  <  9i  <  9  +  \, 


(5.10) 


5-5 


where  g  is  the  average  of  all  the  GIPs  for  the  secondary  data  set  and  is  chosen  large  enough 
such  that  only  20  secondary  data  vectors  satisfy  the  relationship.  A  good  selection  scheme  might 
be  to  rank  in  ascending  order  the  secondary  data  vectors  according  to  the  absolute  value  of  their 
associated  GIP  with  the  mean  GIP  removed.  The  selection  of  the  most  homogeneous  is  then  just 
a  matter  of  choosing  the  first  20  vectors. 

To  understand  why  the  GIP  works  as  a  non-homogeneity  detector,  first  rewrite  Eqn.  (5.9)  in 
the  form 


xf  ^Xi  =  (R»  ^Xi  = 


-  _1 


(5.11) 


where  the  square  root  of  the  matrix  can  be  determined  through  its  eigenvalue/eigenvector  decom¬ 
position.  The  final  term  shown  is  the  output  power  of  a  whitening  filter  operating  on  the  data 
vector.  If  we  denote  the  output  of  the  whitening  filter  as  Xi^  then  the  covariance  matrix  of  this 
output  is 


Ri  =  ^{Xixf}  =  ^{R: 


.  “Xixf  ’ 


}  =  r: 


^  R^Ra 


(5.12) 


where  R^  represents  the  true  covariance  matrix  of  the  data  vector  under  test,  Xi*  Examining  this 
equation  shows  that  as  the  covariance  matrix  estimate  R^  approaches  the  true  covariance  matrix 
Rf,  then 


R.j  ^  Inm  !Rj}  y  R.^.  (5.13) 

This  relationship  further  solidifies  the  whitening  filter  concept.  If  the  data  is  actually  correlated 
according  to  the  estimate  R^,  then  the  data  is  successfully  whitened  to  a  zero  mean  unit  variance 
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random  vector.  Hence,  the  expected  value  of  the  GIF  becomes 


£{9i}  =  £{xfK^Xi}  (5-14) 

where  is  the  whitened  data.  Since  this  data  is  now  white  and  the  data  samples  are  assumed 
independent,  then  the  inner  product  shown  is  the  sum  of  the  unit  variances.  Hence, 

£{9i}  =  £  {lXi[0]|"  +  IXiWl"  +  ■ . .  +  \UNM  -  if  }  =  NM.  (5.15) 

Any  deviation  from  this  value  flags  a  non-homogeneity  in  the  data. 

Using  this  new  knowledge  about  the  homogeneity  of  the  data,  the  MLE  estimator  for  the 
range  gate  can  be  rewritten  as 

Rj  =  j  E  (5.16) 

9ih^ 

where  $  is  the  set  of  J  vectors  satisfying  Eqn.  (5.10).  Although  the  above  formula  only  stipulates 
i  ^  typically  two  guard  cells  on  each  side  of  the  range  cell  of  interest  are  also  removed  from  the 
secondary  data  used  in  the  interference  covariance  matrix  calculation. 

5,4  MCARM  Example 

To  illustrate  the  effectiveness  of  the  non-homogeneity  detector,  an  example  is  offered  from 
the  MCARM  database.  See  the  MCARM  reports  [18]  for  a  description  of  the  database  and  the 
flight  tests  used  to  generate  the  data. 

The  data  file  used  for  this  example  is  rd050465.  Because  of  the  limited  size  of  the  Coherent 
Processing  Interval  (CPI),  the  Factored  Time-Space  (FTS)  algorithm  was  used.  The  number  of 
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range  cells  in  the  CPI  supports  the  sample  size  requirement  for  the  interference  covariance  matrix 
estimate  R. 

A  Blackman  Harris  window  [11]  was  used  for  the  Doppler  filter.  This  window  posesses  excel¬ 
lent  sidelobe  level  characteristics.  The  first  sidelobe  level  is  -92  dB  with  -6  dB/Octave  roll-off  from 
the  first  sidelobe  [11]. 

The  output  of  the  Doppler  filter  is  two  dimensional.  There  is  a  Doppler  filter  bank  for  every 
range  cell.  Within  each  bank  there  are  M  bins,  equal  to  the  number  of  pulses  in  the  CPI.  Although, 
the  output  of  the  filter  could  be  examined  with  a  three  dimensional  plot,  many  of  the  important 
characterstics  are  not  easily  visible.  For  this  reason,  two  figures  are  presented  to  describe  the  output 
of  the  Doppler  filter. 

The  Doppler  filtering  operation  is  performed  for  each  antenna  element.  Therefore,  there  are 
N  Doppler  filters  each  with  a  two  dimensional  output.  If  we  wanted  to  form  a  beam  perpendicular 
to  the  antenna  array  (boresight),  we  could  simply  add  the  output  of  the  Doppler  filters  together. 
However,  the  transmit  direction  for  this  data  file  is  not  boresight,  it  is  19.155®  (reference  the 
MCARM  data  manuals  for  the  coordinate  system  [18]).  Hence,  to  examine  the  output  of  the 
Doppler  filter  at  the  target  location  we  must  form  a  beam  in  the  direction  of  the  target.  The 
target  location  was  conveniently  chosen  to  correspond  to  the  transmit  direction  of  the  radar.  The 
beamforming  operation  is  exactly  the  same  as  that  described  in  Chapter  II.  There  is  no  spatial 
windowing  used  in  the  non-adaptive  beamforming  of  the  Doppler  filters. 

Figure  5.1  shows  the  output  of  the  Doppler  filter  bank  corresponding  to  range  cell  450.  The 
dimension  is  across  the  M  bins  of  the  bank  (there  are  128  pulses  in  the  CPI  and  hence  128  bins  in 
the  Doppler  filter  bank).  The  filtering  operation  was  performed  on  both  the  original  data  file  and 
the  data  file  with  an  artificial  target  inserted.  The  22  (there  are  22  antenna  channels  in  the  CPI) 
Doppler  filters  were  beamformed  to  the  transmit  and  target  azimuth  of  19.155®  with  no  spatial 
windowing  used.  Table  5.1  presents  the  parameters  of  the  target.  The  parameter  of  interest  for 
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Beamformed  Blackman  Harris  Doppler  Filter  Output 


Doppler  Bin 

Figure  5.1  The  Blackman  Harris  Doppler  filter  output  beamformed  to  the  trans¬ 
mit  angle  of  19.155®  for  range  cell  450.  The  artificial  target  is  located  at  Doppler  bin 
65.  The  effect  of  inserting  the  target  into  the  data  is  clearly  visible  by  comparing 
the  two  curves. 
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Table  5.1  Target  parameters  used  for 
the  GIP  non-homogeneity  detector  exam¬ 
ple.  _ 


Variable 

Value 

Range  Cell 

450 

Doppler  Frequency 

1000  Hz 

Doppler  Bin 

65 

Azimuth 

19.155° 

this  plot  is  the  Doppler  bin  location  of  the  target,  bin  65.  Examining  the  figure  shows  the  effect 
of  inserting  the  target  into  the  data,  a  small  spike  at  bin  65.  Without  the  original  data  curve  as  a 
reference,  it  is  impossible  to  determine  the  target  exists  at  this  point. 

The  second  dimension  of  the  Doppler  filter  output  is  across  the  L  (the  number  of  range  cells, 
630)  filter  banks.  We  are  only  interested  in  examining  the  output  of  bank  450,  corresponding  to 
the  location  of  the  target  at  range  cell  450.  Rather  than  plot  across  all  630  range  cells,  a  small 
subset  is  chosen.  In  this  case  the  subset  is  from  cells  420  to  490.  Because  the  target  is  located  at 
Doppler  bin  65,  this  is  the  bin  chosen  for  the  plot  in  Fig.  5.2.  Both  the  original  data  and  the  data 
with  the  artifical  target  is  plotted.  Without  the  original  data  as  a  reference,  detection  of  the  target 
is  not  possible  at  this  point. 

The  amplitude  of  the  target  was  chosen  to  highlight  the  advantage  of  an  adaptive  filter¬ 
ing/beamforming  operation.  Up  to  this  point  in  the  example,  the  presence  of  the  target  cannot 
be  detected.  Typically,  a  detection  scheme  would  use  Fig.  5.2  as  a  method  to  determine  target 
presence.  A  threshold  is  set  and  target  presence  declared  if  the  output  amplitude  of  the  Doppler 
bin  of  interest  at  the  range  cell  of  interest  exceeds  this  threshold.  At  this  point  in  the  example, 
we  have  chosen  an  exceptional  Doppler  filter  window  (the  Blackman  Harris)  and  beamformed  in 
the  direction  of  the  target.  Yet  we  still  cannot  verify  its  presence.  The  only  option  left  is  window¬ 
ing  in  the  spatial  domain  (in  the  non-adaptive  beamforming  operation).  However,  for  comparison 
purposes  to  the  adaptive  FTS  algorithm,  spatial  windowing  will  not  be  used  in  the  example. 
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Doppler  Bin  65  Magnitude  (dB) 


Beamformed  Blackman  Harris  Doppler  Filter  Output 


420  430  440  450  460  470  480  490 

Range  Cell 


Figure  5.2  The  Blackman  Harris  Doppler  filter  output  beamformed  to  the  trans¬ 
mit  azimuth  of  19.155°  for  Doppler  bin  65.  The  artificial  target  inserted  at  range 
cell  450  is  not  threshold  detectable  at  this  point,  there  are  other  peaks  of  higher 
amplitude.  However,  given  the  original  data  as  a  reference  (the  solid  line),  the 
target  is  clearly  visible. 
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FTS  Output 


Figure  5.3  The  FTS  (with  Blackman  Harris  Doppler  filter)  algorithm  output. 
The  target  inserted  at  range  cell  450,  Doppler  bin  65  is  distinctly  shown.  The  AMF 
CFAR  test  statistic,  shown  by  the  dashed  line,  illustrates  a  distinct  improvement 
over  the  raw  FTS  output.  This  statistic  also  allows  the  use  of  a  constant  threshold 
based  on  probability  of  false  alarm. 
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AMF  CFAR  of  ^65  (dB) 


Now  let’s  implement  the  FTS  algorithm  and  examine  the  output  for  the  target.  See  Section  3.5 
for  a  discussion  of  FTS.  No  spatial  windowing  is  used  in  the  implementation  of  the  algorithm  to 
allow  easy  comparison  to  the  non-adaptive  plots  just  presented.  The  Doppler  filter  is  the  same  as 
that  used  above,  the  Blackman  Harris.  The  number  of  secondary  data  vectors  required  (to  give 
an  output  within  3  dB  of  optimum  [15])  for  the  covariance  matrix  estimate  is  twice  the  number 
of  DOF.  Since  the  FTS  algorithm  is  adaptive  in  the  spatial  domain  only,  there  are  22  DOF  (11 
over  11  channel  array)  and  44  secondary  data  vectors  are  required.  For  the  FTS  output  shown 
in  Fig.  5.3,  these  44  vectors  are  chosen  symmetrically  about  the  range  cell  under  test  (using  the 
same  Doppler  bin  from  each  range  cell),  neglecting  two  guard  cells  on  each  side.  For  example,  the 
covariance  matrix  estimate  for  range  cell  450  and  Doppler  bin  65  (the  target  cell  and  bin)  will  use 
the  secondary  data  vectors  from  Doppler  bin  65  and  range  cells  426-447,  453-474. 

As  Fig.  5.3  shows,  the  target  is  clearly  visible  at  range  cell  450.  In  fact,  it  is  approximately  5 
dB  larger  than  the  next  highest  peak.  Threshold  detection  would  work  well  here,  resulting  in  target 
detection.  The  problem  with  the  raw  FTS  output  is  it  does  not  posess  a  Constant  False  Alarm 
Rate  (CFAR)  property.  For  this  reason,  the  output  of  the  algorithm  depends  on  the  noise  and  a 
constant  threshold  cannot  be  set.  Therefore,  several  test  statistics  have  been  introduced  posessing 
the  desirable  CFAR  property.  This  property  allows  a  single  threshold  once  probability  of  false 
alarm  has  been  chosen.  See  Appendix  C  for  a  discussion  of  these  test  statistics  and  more  details 
on  the  theory  behind  them.  For  this  example,  the  Adaptive  Matched  Filter  CFAR  (AMF-CFAR) 
test  was  chosen.  It  is  also  plotted  in  Fig.  5.3,  it’s  amplitude  corresponds  to  the  scale  on  the  right 
side  of  the  graph.  The  embedded  CFAR  characteristic  also  results  in  a  greater  separation  between 
the  target  peak  and  the  next  lower  peak.  It  is  increased  from  approximately  5  dB  to  greater  than 
10  dB. 

Now  let  us  examine  the  output  of  the  FTS  algorithm  when  the  GIF  non-homogeneity  detector 
is  used.  First,  the  interference  covariance  matrix  is  determined  from  the  full  set  of  secondary  data. 
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Secondary  Data  Vector  Sorting 


Figure  5.4  The  GIF  and  sorting  statistic  used  to  order  the  secondary  data  vectors 
according  to  homogeneity.  The  most  homogeneous  vectors  are  closest  to  the  mean 
(the  horizontal  line).  The  sorting  statistic  is  shown  in  Eqn.  (5.17).  This  plot 
illustrates  the  sorting  statistic  results  in  the  correct  ordering  of  the  secondary  data 
vectors  by  GIF  magnitude. 
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In  this  case,  the  secondary  data  set  size  was  set  to  four  times  the  required  number  of  secondary 
data  vectors.  Therefore,  the  secondary  data  set  size  was  set  to  176.  These  176  data  vectors  are 
chosen  symmetrically  about  the  range  cell  under  test.  Figure  5.4  shows  the  ordering  performed 
using  the  GIF.  The  test  statistic  used  to  order  the  vectors  was 

sorting  statistic  =  \gi  -  g\ ,  (5.17) 

where  gi  is  the  GIF  for  the  secondary  data  vector  and  g  is  the  mean  of  all  the  GIFs  within  the 
secondary  data  set.  The  figure  shows  both  the  test  statistic  and  the  true  magnitude  of  the  GIF.  The 
test  statistic  does  result  in  the  correct  ordering  of  the  data  vectors,  since  the  most  homogeneous 
are  clustered  about  the  mean.  Choosing  the  44  secondary  data  vectors  according  to  this  ordering 
results  in  those  closest  to  the  mean  being  chosen  first  as  desired.  Selection  is  merely  a  matter  of 
choosing  the  first  44  secondary  data  vectors. 

Figure  5.5  is  an  alternate  way  of  viewing  the  sorting  process.  The  x-axis  of  the  plot  shows 
the  sorted  data  vector  number.  This  number  is  the  order  in  which  the  data  vectors  are  used  in  the 
MLE  for  the  interference  covariance  matrix  estimate.  The  y-axis  is  the  corresponding  original  data 
vector  number  in  relation  to  the  range  cell  under  test.  The  value  of  0  corresponds  to  the  range  cell 
under  test,  hence  this  data  vector  is  not  used  in  the  MLE  as  illustrated  by  the  absence  of  an  asterisk 
along  this  line  in  the  plot.  The  data  in  the  plot  is  for  range  cell  490.  Similarly,  the  two  guard  cells 
on  each  side  of  the  range  cell  under  test  are  also  not  used  in  the  MLE.  As  the  figure  shows,  the 
use  of  GIF  ordering  results  in  a  drastically  different  set  of  secondary  data  vectors,  since  the  first 
44  of  the  ordered  set  are  chosen  to  generate  the  new  interference  covariance  matrix  R  estimate. 
Also,  the  target  in  range  cell  450  represents  the  most  severe  non-homogeneity  in  the  data,  hence  it 
is  chosen  last  as  indicated  by  the  arrow  in  Fig.  5.5. 

Figure  5.6  shows  the  output  of  the  FTS  algorithm  when  using  the  GIF  non-homogeneity 
detector.  The  target  at  range  cell  450  is  clearly  visible  and  the  difference  between  the  target 
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Original  Data  Vector  Number 


Impact  of  Sorting 


Ordered  Secondary  Data  Vector  Number 


Figure  5.5  A  different  viewpoint  of  the  ordering  process.  The  sorting  data  is  from 
the  interference  covariance  matrix  estimate  for  range  cell  490.  Cell  490  corresponds 
to  0  on  the  y-axis.  This  figure  illustrates  the  impact  of  sorting  the  data  vectors 
according  to  homogeneity.  As  expected,  the  target  represents  the  most  severe  non¬ 
homogeneity  and  hence  it  is  chosen  last  (data  point  shown  by  the  arrow  corresponds 
to  cell  490-40=450.). 


5-16 


Figure  5.6  The  output  of  the  FTS  algorithm  for  Doppler  bin  65  when  using  the 
GIF  non-homogeneity  detector.  The  target  is  clearly  shown  at  range  cell  450  with 
a  difference  of  7  dB  to  the  next  highest  peak.  The  AMF  CFAR  statistic  increases 
this  difference  to  13  dB. 


peak  and  noise  floor  is  considerably  larger  than  without  non-homogeneity  detection.  Figure  5.6 
also  shows  the  output  of  the  FTS  algorithm  when  using  the  AMF  CFAR  statistic.  Even  more 
separation  is  obtained  between  the  target  peak  and  the  noise  floor. 

5.5  Summary 

This  chapter  introduced  the  underlying  architecture  behind  determining  the  interference  co- 
variance  matrix  estimate  R  from  real  world  data.  The  estimate  is  obtained  through  the  use  of  the 
MLE.  However,  the  primary  limitation  to  the  MLE  is  the  assumption  of  iid  data.  This  assumption 
is  not  generally  valid  in  the  real  world. 

Melvin’s  GIF  was  introduced  as  a  method  to  detect  non-homogeneities  in  the  secondary  data. 
The  philosophy  behind  its  operation  was  discussed  in  detail.  As  proof  of  its  operation,  an  example 
using  the  MCARM  database  is  shown  in  Section  5.4.  The  flgures  shown  in  this  section  illustrate 
the  improvement  obtained  by  using  this  non-homogeneity  detector. 

Although  the  example  presented  in  this  chapter  shows  considerable  improvement,  an  even 
greater  improvement  is  obtained  if  a  large  second  target  is  placed  in  the  data.  This  second  target 
corrupts  the  interference  covariance  matrix  estimate  for  the  range  cell  of  the  first,  making  the  orig¬ 
inal  target  undetectable  even  using  the  AMF  CFAR  test  statistic.  Using  the  GIF  non-homogeneity 
detector  results  in  a  target  detection. 

The  most  significant  drawback  to  the  use  of  the  GIF  as  a  non-homogeneity  detector  is  the 
computational  cost  associated.  The  interference  covariance  matrix  must  first  be  estimated  for  the 
entire  secondary  data  set.  Upon  completion  of  this  estimate,  all  of  the  secondary  data  vectors  must 
be  rank  ordered.  The  final  version  of  R  is  obtained  using  only  the  most  homogeneous  members  of 
the  secondary  data  set.  The  computational  cost  associated  with  this  is  immense.  It  severely  offsets 
the  gains  associated  with  implementing  the  reduced  dimension  algorithm. 
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VI.  Analysis  of  the  CSM  Algorithm 


6  A  Introduction 

Previous  chapters  have  explored  the  concepts  of  Space-Time  Adaptive  Processing  (STAP) 
and  extended  the  material  through  many  examples  and  Matlab®  simulations.  The  thesis  now 
focuses  on  satisfying  the  initial  goals  ,  determination  of  the  sample  support  required  for  the  Cross 
Spectral  Metric  (CSM)  algorithm  and  its  susceptibility  to  non-homogeneities  in  the  secondary  data. 

A  second  issue  addressed  is  the  question  of  sample  support  for  the  Maximum  Likelihood  Esti¬ 
mate  (MLE)  of  the  interference  plus  noise  covariance  matrix,  R.  As  discussed  briefly  in  Chapter  V, 
the  analysis  done  by  Reed  [15]  applies  to  STAP  methods  using  Sample  Matrix  Inversion  (SMI), 
There  are  distinct  differences  between  the  SMI  techniques  and  the  CSM  algorithm  that  bring  into 
question  the  application  of  Reed’s  analysis. 

This  chapter  is  organized  into  three  main  sections.  The  first  section  analyzes  the  performance 
of  the  CSM  algorithm  in  terms  of  the  Signal-to-Interference-plus  Noise  Ratio  (SINR),  The  analysis 
is  done  for  a  known  interference  plus  noise  covariance  matrix,  R.  This  analysis  gives  optimum 
performance  curves,  where  optimum  signifies  the  curves  represent  an  SINR  upper  bound  for  the  R 
of  interest. 

These  curves  serve  as  a  foundation  for  comparing  the  results  of  the  next  section.  This  second 
section  attempts  to  define  the  secondary  data  support  required  in  the  MLE  to  obtain  a  certain 
threshold  of  the  optimum  SINR  performance.  The  analysis  is  done  through  the  use  of  Monte  Carlo 
simulation. 

The  third  section  of  this  chapter  examines  the  effects  of  non-homogeneities  in  the  secondary 
data.  As  shown  in  Chapter  V,  the  Factored  Time  Space  (FTS)  algorithm  is  susceptible  to  non¬ 
homogeneities  in  the  secondary  data.  Here,  the  susceptibility  of  the  CSM  algorithm  is  examined 
and  the  Generalized  Inner  Product  (GIP)  non-homogeneity  detector  [14]  is  implemented  with  this 
algorithm. 
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ff.S  SINR  Comparisons  for  Known  Covariance 


6.2,1  Fully  Adaptive.  Let  us  begin  the  analysis  of  the  CSM  algorithm  through  an  ex¬ 
amination  of  SINR  with  known  interference  plus  noise  covariance  R.  The  direct  form  and  FTS 
algorithms  of  Chapter  III  are  used  for  comparison. 

Prom  Eqn.  (3.91),  the  SINR  in  terms  of  the  full  dimensioned  weight  vector  and  R  is 


SINR  = 


|w"v(|^ 
w^Rw  ’ 


(6.1) 


where  is  the  noise  power  per  element,  is  the  target  SNR  per  element,  and  vt  is  the  un¬ 
normalized  space-time  steering  vector.  This  equation  can  be  applied  to  the  CSM  algorithm  if  the 
weight  vector  describing  the  entire  model  is  developed. 

Deriving  a  weight  vector  that  fully  describes  the  CSM  algorithm  is  rather  simple.  We  will 
begin  with  the  full  dimension  model.  Begin  by  writing  the  output  of  the  full  dimensional  CSM 
algorithm  as 


y-d-yh,  (6.2) 

where  y  is  the  output  of  the  CSM  algorithm,  d  is  the  output  of  the  upper  branch,  and  yh  is  the 
output  of  the  lower  branch.  Figure  4.3  depicts  the  model  we  are  currently  working  with.  Prom  the 
development  in  Chapter  IV  and  the  block  diagram  in  Fig.  4.3,  we  can  write  Eqn.  (6.2)  as 

y  = 

=  (8^-w^B)x,  (6.3) 
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where  Wb  is  the  MN  —  1  dimensional  weight  vector  and  s  is  the  normalized  space-time  steering 
vector,  both  from  Section  4.3.  The  equation  now  fits  the  form  of  the  standard  fully  adaptive 
algorithm  y  =  with  the  hermitian  of  the  full  dimensioned  weight  vector  as 

—  Wb  B.  (6.4) 

Substituting  Eqn.  (4.16)  for  Wb  and  taking  the  hermitian  of  both  sides  results  in  the  following  form 
for  the  full  dimension  weight  vector, 


w  =  s  -  B^Rb  ^rbd. 


(6.5) 


We  can  now  substitute  this  weight  vector  description  of  the  fully  adaptive  CSM  algorithm 
into  the  SINK  expression  of  Eqn.  (6.1).  This  substitution  results  in  SINRcsm  full  given  by 


SINRcsm  full  = 


(s-B^Rb  ^rbd)^Vi 


(s  -  BffRb-^rbd)"R  (s  -  B^Rb-^Pbd)  ’ 


(6.6) 


where  the  subscripted  SINRcsm  full  is  to  distinguish  this  expression  from  future  expressions. 

Figure  6.1  shows  a  plot  of  Eqn.  (6.6)  versus  normalized  Doppler  Q,  For  comparison  to  the 
SINK  plots  shown  in  Chapter  III,  both  and  have  been  set  to  1  and  the  interference  plus  noise 
covariance  matrix  R  was  modeled  using  the  same  scenario  as  the  examples  in  Chapter  III.  This 
scenario  includes  two  noise  jammers  (located  at  —40®  and  25®),  clutter,  and  noise.  The  antenna 
array  is  composed  of  18  elements  and  18  pulses  were  used  in  the  CPI. 

It  is  important  to  point  out  some  potential  pitfalls  in  using  Eqn.  (6.6).  When  plotting  the 
equation,  the  steering  vectors  s  (normalized)  and  v  (un-normalized)  are  not  the  only  parameters 
changing  with  normalized  Doppler  Q  (and/or  spatial  frequency  if  a  plot  of  angle  is  desired  also). 
Because  the  blocking  matrix  B  is  defined  such  that  Bs  =  0,  then  every  time  s  changes  the  blocking 
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Fully  Adaptive  CSM  SINK 


Figure  6.1  A  plot  of  the  fully  adaptive  SINK  expression  of  the  CSM  algorithm 
versus  normalized  Doppler.  The  interference  plus  noise  covariance  matrix  R  is  the 
same  as  the  one  used  to  generate  the  SINK  plots  in  Section  3.5  (the  CPI  includes  18 
pulses  and  18  antenna  elements  with  R  including  effects  from  two  jammers,  noise, 
and  clutter).  As  expected,  this  plot  of  the  CSM  algorithm  fully  adaptive  SINK 
matches  the  model  from  Chapter  III. 


matrix  changes  also.  Therefore,  it  is  also  dependent  on  Q,  Furthermore,  the  interference  plus  noise 
covariance  matrix  projection  Rb  is  equal  to  BRB  and  cross  correlation  vector  r^d  is  equal  to  BRs. 
Hence,  Rb  and  rbd  are  also  dependent  on  Q.  This  makes  the  calculation  of  SINRcsm  full  very 
intensive  as  the  inverse  of  Rb  must  be  recalculated  at  each  a). 

All  of  these  pitfalls  have  been  avoided  in  Fig.  6.1.  As  the  figure  depicts,  the  CSM  fully 
adaptive  model  produces  results  very  similar  to  the  model  developed  in  Chapter  III.  The  only 
difference  is  the  peak  SINK  of  25.01  dB,  where  the  previous  model  generated  a  peak  of  25.1  dB. 
The  difference  is  likely  due  to  rounding  errors  incurred  in  the  more  complex  SINK  expression  for 
this  model. 
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The  same  SINK  examination  can  be  computed  for  the  full  dimension  CSM  algorithm  with 
eigenvector  projection  matrix  U  shown  in  Fig.  4.4.  Examination  of  the  figure  shows  the  output  of 


the  model  can  be  written  in  the  same  manner  as  the  above  work, 

y  =  d-yp.  (6.7) 

The  upper  branch  of  the  model  is  exactly  the  same  as  before,  the  projection  of  the  incoming  data 
vector  on  the  normalized  steering  vector  to  give 

y  =  -  Vp-  (6-8) 

The  lower  branch  output  is  the  inner  product  of  p  and  the  weight  vector, 

y  =  s^x  -  Wp^p.  (6.9) 

Working  backwards  through  the  block  diagram  builds  the  expression  allowing  us  to  pull  out  the 
full  dimension  weight  vector, 

y  =  s^x-Wp^U^b 

=  s"x-Wp^U^Bx.  (6.10) 

The  expression  for  the  weight  vector  Wp  was  previously  introduced  in  Eqn.  (4.21), 

Wp  =  U^Wfa.  (6.11) 
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Substituting  this  expression  results  in 


j/  =  s^X-w^UU"Bx.  (6.12) 

Since  the  eigenvector  projection  matrix  U  is  unitary,  the  product  UU^  is  the  identity  matrix. 
Furthermore,  substituting  for  the  weight  vector  Wb  produces 

y  =  (6.13) 

Factoring  out  the  incoming  space-time  snapshot  (or  data  vector)  x  shows  this  expression  is  the 
same  as  that  derived  in  Eqn.  (6.5), 

y=  (s^  -  rgjRb'^B)  x  (6.14) 

which  implies  that 


w  =  s  -  B^Rb  ^rbd. 


(6.15) 


The  two  forms  of  the  full  dimension  CSM  algorithm  have  been  shown  equivalent  in  SINR. 

6,2.2  Partially  Adaptive,  We  have  gained  insight  into  the  CSM  algorithm  through  a 
comparison  of  its  fully  adaptive  SINR  to  the  fully  adaptive  SINR  from  the  model  developed  in 
Chapter  III.  Let’s  gain  further  insight  by  examining  the  SINR  for  the  reduced  dimension  model. 

Referencing  Fig.  4.5  and  Eqn.  (4.30),  the  output  of  the  reduced  dimension  CSM  algorithm 
can  be  written  as 

yr  =  d-yt  =  d-  w^U^h.  (6.16) 
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Substituting  for  d  and  b  in  the  above  equation  and  factoring  out  the  space-time  snapshot  x  results 
in  the  form  needed  to  formulate  a  full  dimension  weight  vector  describing  the  entire  CSM  algorithm, 


(s"-w,^W"b)x.  (6.17) 

Prom  this  equation,  the  full  dimension  weight  vector  for  the  partially  adaptive  CSM  algorithm  is 

w  =  s  -  (6.18) 


Furthermore,  we  can  substitute  Eqn.  (4.26)  to  obtain 


w  =  s  -  B^UAQ^U^Tbd, 


(6.19) 


where  Aq  is  the  diagonal  matrix  of  eigenvalues  corresponding  to  the  O  eigenvectors  used  to  compose 
U.  These  eigenvectors  are  chosen  through  two  popular  methods,  the  cross  spectral  metric  as 
developed  in  Chapter  IV  and  the  principal  components  method.  The  principal  components  method 
simply  chooses  the  eigenvectors  according  to  the  largest  eigenvalue  magnitude.  However,  as  shown 
in  Chapter  IV  this  method  does  not  result  in  the  optimum  SINK. 

The  next  step  is  to  substitute  the  full  dimension  weight  vector  for  the  partially  adaptive  CSM 
algorithm  into  Eqn.  (6.1).  This  substitution  results  in  SINRcsm  part  given  by 


SINRcSM  part  — 


2 

(s-B^WAo'«^rbrf)  v< 

(s  -  B^UAo^U^TbdY  R  (s  -  B^UAo^U^rbd) 


(6.20) 


The  same  pitfalls  exist  here  as  in  the  equation  for  the  fully  adaptive  CSM  algorithm.  The  parameters 
s,  V,  B,  M,  Ao,  and  rbd  all  depend  on  u. 
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SINK  Comparison 


CO 


Figure  6.2  Reduced  dimension  SINK  using  the  cross  spectral  metric  compared  to 
full  dimensionality.  The  dimension  is  18  x  18/16  «  20.  The  full  dimensionality  is 
18  X  18  =  324  DOF. 


Figure  6.2  shows  a  plot  of  Eqn.  (6.20)  versus  normalized  Doppler  a).  When  compared  to 
Fig.  3.20  for  the  FTS  partially  adaptive  algorithm,  the  curve  seems  rather  rough  in  texture.  How¬ 
ever,  this  behavior  does  make  sense.  The  FTS  algorithm  is  adaptive  in  the  spatial  domain  only, 
hence  the  ripples  in  the  curves  shown  in  Fig.  3.20  versus  Q  represent  the  bins  in  the  Doppler  filter. 
The  partially  adaptive  CSM  model  is  adaptive  in  both  the  spatial  and  temporal  domains.  Hence, 
the  DOF  are  chosen  from  both.  This  selection  does  not  necessarily  result  in  a  smooth  pattern  when 
plotted  versus  a).  However,  the  partially  adaptive  CSM  algorithm  chooses  the  DOF  according  to 
the  cross  spectral  metric,  resulting  in  the  optimum  SINR  as  shown  in  Chapter  IV. 

At  this  point,  it  is  interesting  to  compare  the  results  to  the  those  obtained  in  Chapter  III 
for  the  FTS  algorithm.  Figure  6.2  corresponds  to  a  partially  adaptive  algorithm  with  20  DOF. 
Figure  3.20  corresponds  to  the  FTS  algorithm  with  18  DOF.  Comparing  the  two  figures  shows  the 
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effect  of  a  high  quality  Doppler  filter  on  the  FTS  algorithm.  The  CSM  algorithm  is  outperformed 
by  FTS  when  either  the  Blackman  Harris  (-92  dB  sidelobes,  -6  dB/Octave  roll-off)  or  the  Hanning 
(-35  dB  sidelobes,  -18  dB/Octave  roll-off)  windows  are  used  with  the  non-adaptive  Doppler  filter. 
The  FTS  algorithm  performs  better  because  the  CSM  algorithm  simply  cannot  build  the  same 
quality  filter  with  only  20  DOF.  It  must  build  a  two  dimensional  filter  whereas  the  FTS  algorithm 
only  has  to  construct  a  spatial  dimension  filter.  Although  the  CSM  algorithm  doesn’t  produce  the 
highest  peak  SINK,  the  curve  is  more  nearly  horizontal  as  a  function  of  u.  Adding  a  spatial  window 
as  shown  by  Fig.  3.21  does  improve  the  performance  of  the  FTS  algorithm  in  this  regard,  however 
the  maximum  SINK  is  correspondingly  reduced.  Using  the  spatial  tapering,  the  Blackman  Harris 
FTS  algorithm  is  the  only  one  to  outperform  the  CSM  algorithm. 

6,3  Secondary  Data  Support  for  Covariance  Estimation 

As  mentioned  in  Section  5.2,  there  is  some  question  as  to  the  validity  of  applying  Reed’s 
analysis  of  secondary  data  support  to  the  CSM  algorithm.  Reed  has  shown  the  number  of  secondary 
data  vectors  required  to  be  within  3  dB  of  the  optimal  SINR  when  using  the  Sample  Matrix  Inversion 
(SMI)  technique  is  approximately  twice  the  DOF  or  J  =  2iVboF* 

The  analysis  done  by  Reed  is  beyond  the  scope  of  this  thesis,  therefore  secondary  data  support 
requirements  were  determined  experimentally.  The  first  step  in  the  analysis  required  the  production 
of  a  known  covariance  matrix  R.  For  continuity  in  the  thesis,  the  covariance  matrix  used  in  this 
analysis  is  the  same  as  the  one  generated  in  Chapter  III  and  used  in  all  of  the  examples,  with  the 
exception  of  the  MCARM  (real  aircraft  data)  examples  found  in  Chapter  V.  This  interference  plus 
noise  covariance  matrix  contains  the  effects  of  two  jammers  (located  at  —40°  and  25°),  clutter,  and 
noise. 

The  data  is  generated  using  the  standard  technique  for  building  colored  data  from  zero  mean 
unit  variance  Gaussian  samples.  The  first  step  is  generation  of  the  base  samples  with  the  desired 
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Gaussian  distribution. 


X'^A/'(0, 1),  (6.21) 

where  the  notation  shown  is  the  standard  for  a  zero  mean  unit  variance  Gaussian  distribution  and 
X  is  a  column  vector  of  length  MN  (the  length  of  the  space-time  snapshot). 

If  the  desired  covariance  matrix  (in  our  case,  the  one  used  in  Chapter  III)  is  denoted  R,  the 
Cholesky  decomposition  Q  can  be  found  such  that 


R  =  Q^Q. 


(6.22) 


Setting  the  space-time  snapshot  according  to 

X  =  Q^X,  (6.23) 

results  in  the  data  being  colored  according  to  R  as  desired.  The  tilde  above  x  denotes  it  is  diiferent 
than  the  real  space-time  snapshot  x*  Proof  that  the  covariance  matrix  is  indeed  R  is  easily  shown. 
If  the  covariance  matrix  of  the  generated  samples  is  denoted  R,  then 

R  =  £{xx^} 

=  £:{Q^xx^Q} 

=  Q"f{xx^}Q 
=  Q"IQ 

=  Q"Q 

=  R.  (6.24) 


6-10 


Prom  Chapter  V  and  basic  probability  theory,  we  know  the  expected  value  in  Eqn.  (6.24) 
is  simply  the  mean  of  an  infinite  number  of  samples.  This  gave  rise  to  the  Maximum  Likelihood 
Estimator  (MLE)  for  R, 


(6.25) 


where  J  is  the  number  of  secondary  data  vectors  used  in  the  estimate.  As  J  approaches  infinity,  the 
covariance  matrix  estimate  R  will  approach  R.  Using  this  relationship,  we  can  perform  a  Monte 
Carlo  simulation  at  each  secondary  data  set  size  J  and  determine  experimentally  the  required 
secondary  data  support  to  achieve  a  certain  SINK. 

As  verification  of  the  Monte  Carlo  simulation  process,  Fig.  6.3  is  oflFered.  The  FTS  algorithm 
is  used  with  artificially  generated  space-time  snapshots.  This  algorithm  is  adaptive  in  the  spatial 
domain  only.  Since  there  are  18  antenna  elements  in  this  simulation,  there  are  18  DOF.  Reed’s  rule 
predicts  approximately  2iSrDOF  secondary  data  vectors  are  required  to  achieve  performance  within 
3  dB  of  optimal  SINK. 

In  the  experiment,  there  were  20  trials  for  each  secondary  data  set  size.  The  solid  line 
depicts  the  experimental  SINK  Loss,  Lsinr,  for  the  secondary  data  set  size  shown  on  the  x-axis 
and  is  the  mean  of  the  values  returned  from  each  of  the  20  trials.  The  averaging  technique  used 
is  the  arithmetic  mean,  signifying  the  mean  of  the  samples  is  calculated  in  the  native  units.  The 
conversion  to  dB  is  done  after  the  mean  has  been  calculated.  The  standard  deviation  of  these  20 
trials  is  depicted  by  the  dashed  line.  The  standard  deviation  is  also  calculated  in  the  native  units 
and  then  converted  to  dB.  The  solid  line  in  the  figure  reaches  3  dB  of  optimal  with  secondary 
data  support  between  30  and  35  vectors.  As  expected,  these  results  coincide  very  well  with  Reed’s 
prediction  of  2  x  18  —  3  =  33. 

Figure  6.4  shows  a  Monte  Carlo  simulation  for  the  full  dimension  CSM  algorithm.  Again, 
the  solid  line  depicts  the  average  Lsinr  for  the  20  trials  and  the  dashed  line  depicts  the  standard 
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mean  Lsinr  (dB) 


FTS  Algorithm  with  18  DOF 


20  40  60  80  100  120  140  160  180  200 


J 

Figure  6.3  The  experimental  Tsinr  curve  versus  secondary  data  set  size  for  the 
FTS  algorithm.  There  are  18  DOF,  hence  Reed  predicts  33  secondary  data  vectors 
to  achieve  performance  within  3  dB  of  optimal  SINK.  The  Monte  Carlo  simulation 
indicates  between  30  and  35  secondary  data  vectors  axe  required,  agreeing  very  well 
with  Reed’s  prediction  of  33.  The  analysis  was  done  for  Doppler  bin  9  (of  18  bins), 
no  windows/tapers  were  used. 
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Std  Dev  (dB) 


deviation.  Before  examining  the  results  shown  in  the  figure,  an  explanation  of  the  third  curve  is 
in  order.  This  curve  is  characterized  by  the  dash  dot  sequence  and  depicts  the  ratio  of  the  mean 
SINK  to  the  standard  deviation  of  these  experiments.  When  using  the  dB  values,  this  ratio  simply 
becomes  the  standard  deviation  subtracted  from  the  mean.  The  ratio  provides  an  indicator  of  the 
quality  of  the  Monte  Carlo  simulation. 

Typically,  a  ratio  of  8  to  10  dB  signifies  enough  trials  were  performed  in  the  simulation.  Higher 
ratios  correspond  to  higher  confidence  in  the  results.  According  to  Fig.  6.4,  the  ratio  indicates  the 
20  trials  used  in  the  Monte  Carlo  simulation  were  adequate  for  this  experiment.  When  the  ratio  is 
smaller  as  in  Fig.  6.7,  the  Tsinr  curve  is  not  as  smooth.  This  smoothness  is  another  indication  of 
the  quality  of  the  experiment.  Large  differences  in  mean  Lsinr  as  the  sample  support  is  increased 
suggest  not  enough  trials  were  done  in  the  experiment.  Although  Fig.  6.7  exhibits  a  small  ratio, 
the  mean  Lsinr  values  do  not  exhibit  exaggerated  differences  as  the  secondary  data  support  is 
increased.  Therefore,  the  number  of  trials  was  not  increased  for  the  experiments.  The  ratio  curve 
was  omitted  from  Fig.  6.3  because  the  values  were  on  the  order  of  +30  dB,  well  above  the  scale 
used  for  the  standard  deviation  curve  making  it  undesirable  to  add  this  third  curve  to  the  plot. 

Continuing  the  analysis  of  Fig.  6.4,  the  optimal  SINK  (known  covariance)  for  full  dimension¬ 
ality  was  shown  in  Fig.  6.1  as  25  dB.  The  marked  point  in  Fig.  6.4  represents  3  dB  within  this 
optimal  SINR  of  25  dB.  As  expected  for  the  full  dimension  model,  this  point  corresponds  to  a  data 
set  size  of  approximately  twice  the  DOF,  or  2  x  324  =  648.  Therefore,  approximately  650  secondary 
data  vectors  are  needed  to  estimate  the  covariance  matrix  and  retain  performance  within  3  dB  of 
the  optimal  SINR.  The  Monte  Carlo  simulation  results  in  the  conclusion  that  Reed’s  rule  holds  for 
the  CSM  algorithm  when  using  full  dimensionality. 

The  next  step  in  the  analysis  is  Monte  Carlo  experiments  for  reduced  dimension  CSM  algo¬ 
rithms.  However,  in  order  to  calculate  Lsinr  we  must  first  determine  the  optimum  SINR  (known 
covariance)  for  the  reduced  dimensions  of  interest.  This  calculation  is  required  because  the  opti- 
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Monte  Carlo  Analysis  for  Pull  Dimensionality 


Figure  6.4  The  Lsinr  curve  versus  secondary  data  set  size  for  the  full  dimensioned 
CSM  model.  As  expected  for  full  dimensionality,  it  is  approximately  twice  the  DOF 
(2  X  324  =  648). 
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Std.  Dev.  (dB),  Ratio  of  Mean  SINK  to  Std.  Dev.  (dB) 


Optimum  SINK  Curve 


Figure  6.5  The  optimum  SINK  curve  (known  covariance)  versus  weight  vector 
dimension  (A^dof)*  The  SINK  curve  illustrates  as  few  as  65  DOF  are  required  to 
reach  the  maximum  optimal  SINK. 


mum  SINK  will  experience  roll-off  as  the  dimension  of  the  problem  is  reduced.  Figure  6.5  presents 
the  optimum  SINR  curve. 

Figure  6.5  is  the  reference  needed  to  determine  the  Lsinr  resulting  from  using  the  interference 
plus  noise  covariance  matrix  estimate.  The  figure  also  gives  insight  into  the  reason  behind  reduced 
dimension  STAP  algorithms.  The  optimum  SINR  curve  rapidly  approaches  the  maximum  SINR  of 
25  dB  with  far  fewer  than  MN  DOF.  In  this  example,  the  same  interference  plus  noise  covariance 
matrix  is  being  used  as  in  all  previous  examples.  This  covariance  matrix  corresponds  to  18  x  18,  or 
324  DOF.  The  Ndof  required  to  reach  the  maximum  optimum  SINR  is  approximately  65,  much 
less  than  the  full  dimension  of  324. 

Figure  6.6  is  a  plot  of  the  experimental  Lsinr  versus  secondary  data  set  size  J  for  the  partially 
adaptive  CSM  algorithm  using  65  DOF.  The  solid  line  indicates  the  mean  Lsinr  for  the  20  trials 
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Monte  Carlo  Analysis  for  65  DOF 


Figure  6.6  The  Lsinr  curve  versus  secondary  data  set  size  for  the  partial  dimen¬ 
sion  CSM  algorithm  using  65  DOF.  The  labeled  point  represents  3  dB  of  optimal 
SINK  (25  dB).  Contrary  to  Reed’s  analysis,  the  number  of  secondary  data  vectors 
required  is  greater  than  twice  the  DOF  (2  x  65  =  130).  This  reduced  dimension 
model  requires  160  secondary  data  vectors  to  achieve  3  dB  of  optimal  performance. 
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Std.  Dev.  (dB),  Ratio  of  Mean  SINK  to  Std.  Dev.  (dB) 


Monte  Carlo  Analysis  for  45  DOF 


J 


Figure  6.7  The  Lsinr  curve  versus  secondary  data  set  size  for  the  partial  dimen¬ 
sion  CSM  algorithm  using  45  DOF.  The  labeled  point  represents  4  dB  of  optimal 
SINK  (20  dB).  Contrary  to  Reed’s  analysis,  the  curve  does  not  reach  3  dB  of  opti¬ 
mal.  This  reduced  dimension  model  requires  420  secondary  data  vectors  or  9  times 
the  DOF  to  achieve  4  dB  of  optimal. 

at  each  secondary  data  set  size.  The  dashed  line  corresponds  to  the  standard  deviation  of  these  20 
trials.  The  data  set  size  was  increased  in  steps  of  10.  As  before,  the  dash-dot  line  corresponds  to 
the  ratio  of  the  mean  SINK  to  the  standard  deviation.  The  ratio  curve  shown  in  the  figure  indicates 
the  trial  size  of  20  was  adequate  for  the  simulation. 

The  results  of  this  Monte  Carlo  simulation  indicate  160  secondary  data  vectors  are  required 
to  achieve  SINR  performance  within  3  dB  of  optimal.  This  performance  difiers  from  Reed’s  rule 
stating  the  secondary  data  support  should  be  approximately  twice  the  DOF.  Using  the  2iVDOF  rule 
implies  only  2  x  65  =  130  secondary  data  vectors  are  needed,  resulting  in  the  loss  of  approximately 
1  dB  in  SINR.  The  Monte  Carlo  simulation  indicates  approximately  2.5  secondary  data  vectors  are 
required  for  operation  within  3  dB  of  optimal. 
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Comparison  of  Experimental  Lsinr  Curves 


Figure  6.8  Overlay  of  the  Lsinr  curves  for  reduced  dimensions  below  the  dimen¬ 
sion  of  the  interference  subspace.  Performance  relative  to  optimal  is  reduced  as  the 
DOF  are  reduced. 


The  next  simulation  reduces  the  DOF  to  45.  Fig.  6.5  shows  the  optimum  SINK  (known 
covariance)  is  now  20  dB.  Therefore,  the  curve  needs  to  converge  to  17  dB.  The  Monte  Carlo 
simulation  in  Fig.  6.7  shows  the  Lsinr  curve  never  reaches  this  point.  To  achieve  within  4  dB  of 
optimal  SINK  requires  secondary  data  support  of  size  420.  This  translates  to  420/45  «  9  times  the 
DOF! 

Overlaying  Monte  Carlo  simulation  results  for  several  different  reduced  dimensionalities  shows 
the  trend  of  interest.  Figure  6,8  is  a  plot  the  experimental  SE^R  Loss  (Lsinr)  curves  for  reduced 
dimension  applications  with  the  DOF  below  the  dimension  of  the  interference  subspace.  Each 
curve  has  been  referenced  to  the  appropriate  optimum  SINK  values  from  Fig.  6.5.  The  dimension 
of  the  interference  subspace  was  determined  through  the  use  of  two  tools,  the  Singular  Value 
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Decomposition  (SVD)  of  the  interference  plus  noise  covariance  matrix  R  and  plots  of  the  Cross 
Spectral  Metric  (CSM). 

A  note  about  the  curves  in  Fig.  6.8  is  in  order.  The  sharp  drop-off  each  curve  exhibits 
corresponds  to  the  point  where  J  is  less  than  iVooF-  This  drop-off  occurs  because  the  number 
of  secondary  data  vectors  used  in  the  interference  plus  noise  covariance  matrix  estimate  R  must 
be  at  least  equal  iVooF-  This  requirement  is  a  result  of  the  CSM  algorithm  using  an  eigenvalue 
decomposition  to  reduce  the  dimensionality  of  the  problem.  If  only  J  secondary  data  vectors  are 
used  to  estimate  R  and  J  <  iVboF?  then  only  the  first  J  eigenvalues  are  nonzero.  Referencing 
Eqn.  (6.19),  we  see  the  inverse  of  Aq  (the  diagonal  matrix  of  eigenvalues  corresponding  to  the 
eigenvectors  chosen  by  the  CSM)  is  required.  If  there  are  only  J  nonzero  eigenvalues,  then  Aq 
contains  eigenvalues  that  are  either  equal  to  zero  or  very  small.  Therefore,  the  matrix  is  very  close 
to  singular  and  cannot  be  inverted,  resulting  in  the  sharp  drop-off  reflected  in  each  Lsinr  curve 
shown  in  the  graphs. 

Figure  6.9  shows  the  magnitude  of  the  singular  values  and  the  CSM.  The  upper  ledge 
represents  the  dimension  of  the  interference  subspace  and  the  lower  ledge  is  the  noise  (uncorre¬ 
lated/white).  This  plot  reinforces  the  fact  that  it  is  the  noise  that  makes  R  invertible.  Without  it, 
the  rank  of  R  is  equal  to  the  dimension  of  the  interference  subspace,  which  is  less  than  MN. 

Figure  6.10  is  simply  an  expanded  view  of  Fig.  6.9.  This  allows  a  more  accurate  estimation 
of  the  dimension  of  the  interference  subspace.  From  the  figure,  the  dimension  is  approximately  70. 

These  two  figures  add  considerable  insight  into  the  underlying  architecture  of  the  CSM  algo¬ 
rithm.  The  model  chooses  the  DOF  that  cancel  the  interference,  not  the  noise,  taking  advantage 
of  the  statistical  properties  of  the  incoming  data.  Combining  the  knowledge  gained  from  Figs.  6.9 
and  6.10  with  Fig.  6.8  results  in  the  prediction  that  ATbof  should  be  equal  to  the  dimension  of 
the  interference  subspace  for  optimum  performance  in  situations  of  unknown  covariance,  although 
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Singular  Values  of  R 


Figure  6.9  Singular  value  and  cross  spectral  metric  magnitudes  of  the  known 
interference  plus  noise  covariance  matrix.  The  singular  values  represent  the  power 
contained  in  the  interference  subspace.  This  figure  illustrates  the  most  power  is 
contained  in  the  clutter  and  jamming  interference  subspace.  The  noise  subspace 
represents  very  little  of  the  overall  power. 
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Cross  Spectral  Metric  Magnitude  (dB) 


Singular  Values  of  R 


Singular  Value  Number/CSM  Number 


Figure  6.10  Closer  examination  of  the  most  significant  singular  values  indicates 
the  most  power  is  contained  in  the  first  70,  a  measure  of  the  dimension  of  the 
interference  subspace.  This  plot  coincides  with  Fig.  6.5  in  illustrating  the  DOF 
required  to  reach  optimum  is  equal  to  the  dimension  of  the  interference  subspace. 
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obviously  determining  the  dimension  of  the  interference  subspace  for  real  time  operation  is  very 
expensive. 

Going  back  to  the  Monte  Carlo  analysis,  Fig.  6.8  also  shows  an  unexpected  trend  in  the  data. 
This  trend  is  in  direct  conflict  to  Reed’s  rule.  As  iVboF  is  reduced,  the  required  secondary  data 
set  size  J  to  achieve  a  constant  measure  of  the  optimum  performance  increases.  For  example, 
performance  within  3  dB  of  optimal  is  obtained  with  160  secondary  data  vectors  for  65  DOF  while 
175  are  required  when  the  A^dof  is  reduced  to  60.  Further  reduction  of  JVdof  results  in  SINK 
performance  that  does  not  reach  3  dB  of  optimal.  Each  experimental  Lsinr  curve  is  referenced  to 
the  optimum  SINK  corresponding  to  ATdof- 

The  figure  also  shows  a  distinct  thresholding  effect  occurring  when  iVboF  is  decreased  to 
values  below  55.  These  curves  are  grouped  together  with  all  producing  performance  within  5  or 
6  dB  of  optimal  using  approximately  the  same  number  of  secondary  data  vectors.  These  results 
contradict  the  primary  goal  of  reducing  the  DOF  to  gain  a  corresponding  decrease  in  secondary 
data  support. 

Figure  6.11  overlays  the  experimental  Lsinr  curves  for  DOF  above  the  dimension  of  the 
interference  subspace.  Again,  the  curves  illustrate  the  sharp  drop-off  in  performance  when  J  is 
decreased  below  iVboF*  Contrary  to  the  overlay  shown  in  Fig.  6.8,  increasing  the  DOF  does  not 
increase  the  performance  of  the  algorithm.  Performance  is  significantly  decreased  when  J  is  equal 
to  MAT,  or  324.  Maximum  performance  is  again  obtained  when  ATdof  is  equal  to  the  dimension  of 
the  interference  subspace.  Note  that  once  J  is  increased  beyond  MAT,  the  curves  converge  to  the 
same  point. 


6.3.1  Hypothesis  Behind  Poor  Performance.  The  Monte  Carlo  simulation  has  shown  the 
performance  of  the  CSM  algorithm  when  estimating  the  covariance  is  maximized  when  ATdof  is 
equal  to  the  dimension  of  the  interference  subspace.  Both  Figs.  6.8  and  6.11  support  this  statement. 
The  poor  performance  of  the  algorithm  when  ATdof  is  not  equal  to  the  dimension  of  the  interference 
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Comparison  of  Experimental  Lsinr  Curves 


Figure  6.11  Overlay  of  the  experimental  SINK  curves  for  reduced  dimensions 
above  the  dimension  of  the  interference  subspace.  Performance  relative  to  optimal 
is  reduced  as  the  DOF  are  increased. 
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subspace  is  likely  due  to  calculation  noise  in  the  eigenvalues  of  R.  Small  secondary  data  set  sizes 
result  in  an  estimated  interference  plus  noise  covariance  matrix  R  of  less  than  full  rank.  In  the 
examples  used,  this  statement  is  true  for  any  secondary  data  set  size  less  than  18  x  18,  or  324. 

As  corroborated  by  Fig.  6.9,  the  first  65  to  70  singular  values  are  much  larger  than  any 
subsequent  singular  values.  The  singular  values  are  directly  related  to  the  eigenvalues  [20],  hence 
Fig.  6.9  gives  appropriate  insight.  The  largest  eigenvalues  are  not  affected  as  greatly  by  the  varying 
secondary  data  set  sizes.  However,  the  smaller  eigenvalues  are  greatly  affected. 

Referencing  Eqn.  (4.71),  the  CSM  is  defined  as 


VXi  ’ 


(6.26) 


where  r^d  is  the  cross  correlation  between  the  vector  b  and  the  desired  scalar  d,  i/*  is  the 
eigenvector,  and  A*  is  the  corresponding  eigenvalue.  The  eigenvalue  error  will  be  small  relative  to 
the  magnitude  of  the  largest  eigenvalues,  however  the  error  is  large  relative  to  the  small  eigenvalues. 

This  error  in  the  small  eigenvalues  results  in  an  inherent  instability  in  the  CSM  algorithm 
when  JVdof  is  increased  to  include  the  DOF  corresponding  to  small  eigenvalues.  The  dimension  of 
the  interference  subspace  marks  the  point  at  which  the  eigenvalues  become  much  smaller.  A  small 
change  in  secondary  data  set  size  can  result  in  a  drastically  different  W,  the  dimension  reducing 
matrix  of  eigenvectors.  Remember,  the  eigenvectors  comprising  U  are  chosen  according  to  the 
CSM. 

In  the  examples  shown,  the  dimension  of  the  interference  subspace  was  approximately  65. 
This  is  exactly  the  point  at  which  the  performance  of  the  CSM  algorithm  diminishes,  as  shown  by 
Fig.  6.11.  Once  the  secondary  data  set  size  J  is  increased  beyond  MN  (324  in  the  examples)  the 
performance  begins  to  improve  again.  This  is  explained  through  an  examination  of  the  reduced 
dimension  weight  vector  shown  in  Eqn.  (6.19). 
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Let’s  continue  the  analysis  by  examining  the  performance  of  the  Lsinr  curves  in  a  manner 
from  left  to  right  or  increasing  J.  When  J  <  MN^  there  are  only  J  nonzero  eigenvalues  because 
R  (of  dimension  MN  x  MN)  is  not  full  rank.  As  long  as  JVboF  <  Ji  Aq  contains  only  large 
eigenvalues.  The  matrix  is  easily  inverted.  However,  increasing  Nbof  beyond  J  forces  Aq  to 
contain  eigenvalues  that  are  very  close  to  zero.  Inverting  this  matrix  is  difficult  and  prone  to  errors. 
This  is  shown  by  the  sharp  drop-off  each  Lsinr  curve  exhibits  when  J  <  iVboF- 

When  iVooF  is  equal  to  the  dimension  of  the  interference  subspace,  the  Lsinr  performance  is 
maximum.  This  is  explained  by  again  looking  at  the  eigenvalues.  If  R  was  full  rank,  the  eigenvalue 
decomposition  would  contain  MN  nonzero  eigenvalues.  However,  since  J  is  still  less  than  MAT, 
the  matrix  is  not  full  rank  and  contains  only  J  nonzero  eigenvalues.  The  rest  are  either  zero  or 
very  close  to  zero  depending  on  the  current  random  realization  of  R.  These  small  eigenvalues 
create  calculation  noise  in  the  CSM.  When  ATdof  =  the  CSM  algorithm  is  not  forced  to  choose 
eigenvectors  lying  in  the  noise  subspace  and  it  gets  to  use  all  of  the  eigenvectors  lying  in  the 
interference  subspace.  Hence,  the  performance  of  the  algorithm  is  the  best. 

Choosing  J  <  MN  and  increasing  iVboF  such  that  it  is  larger  than  the  dimension  of  the 
interference  subspace  forces  the  CSM  algorithm  to  choose  eigenvectors  lying  in  the  noise  subspace. 
These  eigenvectors  are  characterized  by  the  extremely  small  eigenvalues.  These  eigenvalues  then 
become  a  part  of  Aq.  When  this  matrix  is  inverted,  the  small  eigenvalues  result  in  the  interference 
estimate  being  larger  than  it  should  be  and  a  corresponding  increase  in  Lsinr*  The  Lsinr  increases 
as  ATdof  increases  because  more  of  the  small  eigenvalues  are  becoming  part  of  Ao  and  corrupting 
the  interference  estimate.  This  is  shown  by  the  results  in  Fig.  6.11. 

Once  J  has  reached  MAT,  the  interference  plus  noise  covariance  matrix  estimate  R  has  reached 
full  rank  and  none  of  the  eigenvalues  are  zero,  although  errors  still  exist  because  of  the  small 
secondary  data  set  size.  The  Lsinr  curves  begin  the  expected  progression  to  optimal  performance 
as  J  is  increased. 
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Eigenvector  Ordering  for  U 


Figure  6.12  The  GSM  ordering  of  the  eigenvectors  for  a  secondary  data  set  size  J 
equal  to  323,  one  less  than  MN,  The  ordering  is  significantly  different  for  estimated 
covariance,  R.  Changing  J  results  in  different  ordering.  The  eigenvectors  in  the 
lower  left  corner  lie  in  the  interference  subspace. 


There  is  still  some  question  as  to  why  the  performance  of  the  GSM  algorithm  experienced 
decreasing  performance  with  iVboF  below  the  dimension  of  the  interference  subspace.  Figure  6.12 
gives  insight  into  this  problem.  The  figure  represents  the  order  in  which  the  DOF  are  chosen  by 
the  GSM  algorithm  for  two  different  scenarios. 

The  first  scenario  assumes  a  known  interference  plus  noise  covariance  matrix,  R.  These  points 
can  be  considered  the  optimal  ordering.  The  a:-axis  is  the  order  in  which  the  eigenvectors  are  chosen 
and  the  ^-axis  is  the  original  order  they  were  in.  For  example,  the  first  eigenvector  or  DOF  in  the 
optimal  ordering  was  originally  the  second  one.  It^s  location  on  the  plot  is  position  (1,2).  This 
means  the  second  DOF  has  the  greatest  positive  effect  (according  to  the  GSM)  on  the  SINK  than 
any  of  the  other  DOF. 
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The  second  scenario  in  the  plot  is  for  an  estimated  interference  plus  noise  covariance  matrix 
R.  The  data  was  generated  in  the  same  manner  as  in  the  Monte  Carlo  experiments,  using  the 
coloring  process  explained  in  that  section.  Only  323  secondary  data  vectors  were  used  in  this 
particular  estimate. 

As  Fig.  6.12  shows,  the  ordering  when  using  R  is  drastically  diiferent.  The  first  65  eigenvectors 
are  the  same  as  indicated  by  the  clustering  in  the  lower  left  corner  of  the  plot.  However,  the  order 
they  are  chosen  has  changed.  This  cluster  contains  the  eigenvectors  that  lie  in  the  interference 
subspace.  As  long  as  iVboF  is  greater  than  65,  the  ordering  of  these  first  data  vectors  doesn’t 
matter  because  the  algorithm  is  using  all  the  vectors  lying  in  the  interference  subspace.  Decreasing 
A^dof  below  65  places  emphasis  on  the  ordering.  This  explains  the  decreased  performance  when 
ATdof  is  below  the  dimension  of  the  interference  subspace,  shown  by  Fig.  6.8. 

6,J^  Non-Homogeneities  and  the  CSM  Algorithm 

In  Chapter  V,  the  impact  of  non-homogeneities  in  the  secondary  data  on  the  output  of 
the  Factored  Time-Space  (FTS)  algorithm  was  explored  in  detail.  The  Generalized  Inner  Product 
(GIP)  [14]  was  introduced  in  an  effort  to  alleviate  this  impact.  The  GIP  was  shown,  using  MCARM 
data  containing  actual  airborne  radar  measurements,  to  result  in  significant  improvement  in  the 
output  of  the  FTS  algorithm.  This  improvement  consisted  of  increasing  the  separation  between 
the  target  peak  and  largest  noise  level.  The  chapter  also  explored  the  mathematical  architecture 
behind  the  GIP  non-homogeneity  detector. 

In  this  chapter,  an  analysis  of  the  CSM  algorithm  is  being  presented.  The  next  logical  step 
in  this  analysis  is  an  examination  of  the  effects  of  non-homogeneities  in  the  secondary  data  on  the 
CSM  algorithm. 

The  examination  begins  exactly  as  the  Monte  Carlo  simulation  for  required  secondary  data 
support  began.  A  Cholesky  decomposition  of  the  desired  interference  plus  noise  covariance  matrix 
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R  gives  a  method  to  color  artificially  generated  white  data.  This  method  is  explained  in  Section  6.3. 
The  same  R  is  used  as  in  all  previous  examples. 

The  simulation  begins  by  generating  one  random  realization  of  a  CPI  datacube  through  the 
coloring  process  explained  in  Section  6.3.  Using  this  one  realization  for  the  three  different  cases 
of  interest  allows  a  direct  comparison  of  the  results.  The  three  ca^es  of  interest  include  operating 
the  CSM  algorithm  on  data  with  only  a  target,  data  with  a  target  and  one  non-homogeneity, 
and  the  CSM  algorithm  with  GIP  non-homogeneity  detection  on  data  with  a  target  and  one  non¬ 
homogeneity.  The  purpose  of  this  section  is  not  to  generate  a  Monte  Carlo  simulation.  The  analysis 
only  requires  generation  of  one  set  of  data  to  determine  the  effects  of  non-homogeneities. 

With  this  goal  in  mind,  a  single  target  was  placed  in  the  CPI  datacube.  This  target  was 
placed  at  range  cell  (?)  100.  To  simulate  a  non-homogeneity  in  the  secondary  data,  a  second  target 
was  placed  at  range  cell  50.  This  second  target  is  identical  to  the  first,  with  the  exception  of  the 
range  cell  location. 

This  non-homogeneity  is  a  worst  case  scenario  for  the  CSM  algorithm.  It  represents  worst 
case  because  the  non-homogeneity  is  the  same  angle  and  Doppler  as  the  desired  target /steering 
vector.  The  CSM  algorithm  cannot  filter  any  of  it  out.  The  gains  seen  here  are  due  solely  to  the 
use  of  the  GIP  detection  scheme. 

Figure  6.13  illustrates  the  output  of  the  CSM  algorithm  both  with  and  without  the  non¬ 
homogeneity  in  the  data.  The  interference  plus  noise  covariance  matrix  estimate  R  was  found 
using  a  symmetric  window  of  secondary  data  vectors  about  the  range  cell  of  interest.  The  CSM 
algorithm  was  programmed  to  use  45  DOF  and  200  secondary  data  vectors.  The  choice  of  200 
secondary  data  vectors  was  made  by  referencing  Fig.  6.8.  This  choice  of  secondary  data  support 
should  give  performance  within  approximately  5  dB  of  optimal  according  to  the  Monte  Carlo 
simulations.  Of  course,  choosing  65  DOF  results  in  better  performance  (within  «  2  dB  of  optimal 
SINK). 
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Figure  6.13  The  impact  of  non-homogeneities  in  the  secondary  data  on  the  CSM 
algorithm.  This  example  shows  the  results  of  inserting  a  target  in  range  cell  50 
identical  to  desired  target  at  range  cell  100.  The  CSM  algorithm  utilizes  45  DOF 
and  R  is  built  with  200  symmetric  secondary  data  vectors. 


Non-Homogeneity  Detection 


Figure  6.14  The  impact  of  using  the  GIF  non-homogeneity  detector  with  the 
CSM  algorithm.  The  curve  is  the  output  of  the  GSM  algorithm  using  45  DOF  and 
200  secondary  data  vectors. 


The  output  of  the  CSM  algorithm  clearly  shows  the  target  at  range  cell  100.  The  non¬ 
homogeneity  has  not  been  inserted  into  the  data  resulting  in  acceptable  performance  when  using  a 
symmetric  window  of  secondary  data  vectors.  However,  introducing  the  non-homogeneity  degrades 
the  performance  of  the  algorithm.  This  is  shown  by  the  second  curve  in  the  graph.  The  separation 
from  target  peak  to  next  highest  peak  has  been  reduced  by  approximately  3  dB.  Since  this  output 
is  typically  threshold  detected,  this  is  a  significant  decrease. 

Figure  6.14  illustrates  the  positive  impact  the  Gff  non-homogeneity  detector  has  on  the  CSM 
algorithm.  The  separation  between  target  and  next  highest  peak  has  been  increased  from  5  dB  to 
8  dB.  This  is  an  improvement  of  60%.  Simply  out  of  curiosity,  the  experiment  was  repeated  after 
modifying  the  CSM  algorithm  to  make  use  of  65  DOF. 


6-30 


Non-Homogeneity  Detection 


Figure  6.15  The  impact  of  using  the  GIF  Non-Homogeneity  detector  with  the 
GSM  algorithm.  The  curve  is  the  output  of  the  GSM  algorithm  using  65  DOF  and 
200  secondary  data  vectors. 


Figure  6.15  shows  the  results  of  this  second  experiment.  Glearly,  the  GIF  non-homogeneity 
detector  is  a  significant  improvement.  The  separation  between  target  and  next  highest  peak  has 
increased  from  4  dB  to  7  dB.  An  increase  of  75%. 

There  is  one  major  limitation  of  the  Gff.  It  requires  the  inverse  of  a  covariance  matrix 
estimate  based  on  the  full  secondary  data  set,  R^.  To  do  the  inverse,  Rj,  must  be  full  rank.  This 
is  a  limiting  factor  when  using  the  GIF.  The  GSM  algorithm  requires  the  estimation  of  a  full 
dimensioned  R,  MAT  x  MAT.  Hence,  R^  is  also  MN  x  MN.  This  translates  to  a  secondary  data 
support  requirement  of  at  least  MN  vectors  for  implementation  of  the  GIF  detection  scheme.  This 
is  in  direct  confiict  to  the  primary  goal  of  reducing  secondary  data  support.  For  the  experiments  in 
this  section,  the  secondary  data  set  size  was  increased  to  324  data  vectors  to  ensure  the  invertibility 
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of  Rs.  Symmetric  windowing  used  200  secondary  data  vectors  symmetric  about  the  range  cell  of 
interest  while  the  GIP  detection  equipped  algorithm  used  the  200  most  homogeneous  secondary 
data  vectors. 

6.5  Summary 

The  analysis  performed  in  this  chapter  has  uncovered  several  interesting  pieces  of  information 
regarding  the  GSM  algorithm.  Expressions  for  the  SINK  were  developed  for  both  full  and  reduced 
dimension  cases.  Plots  of  these  expressions  using  known  covariance  served  to  gain  insight  into 
the  operation  of  the  GSM  algorithm.  Similar  plots  (Ghapter  III)  for  the  FTS  algorithm  allowed 
comparison  between  the  algorithms. 

As  mentioned  previously,  the  GSM  algorithm  requires  the  estimation  of  a  full  dimensioned 
interference  plus  noise  covariance  matrix.  There  is  no  gain  over  the  original  fully  adaptive  model 
in  this  respect.  Furthermore,  the  eigenvalue/eigenvector  decomposition  of  this  matrix  must  be 
performed.  Obviously,  this  algorithm  does  not  enjoy  the  benefits  of  working  in  smaller  dimensions 
such  as  FTS  does.  Galculation  times  are  much  greater  for  this  algorithm  than  for  FTS.  For  this 
reason,  the  GSM  algorithm  will  likely  remain  a  theoretical  tool. 

Although  the  GSM  algorithm  is  computationally  expensive,  it  does  give  optimum  performance 
when  the  dimension  reducing  matrix  is  composed  of  eigenvectors.  This  was  shown  in  Ghapter  IV. 
The  point  to  remember  is  the  dimension  reducing  matrix  does  not  have  to  be  composed  of  eigen¬ 
vectors.  There  may  be  another  matrix  that  gives  better  SINK  performance. 

The  eigenvalue/eigenvector  decomposition  of  R  in  the  GSM  algorithm  brought  into  question 
the  validity  of  Reed’s  rule.  In  an  effort  to  answer  this  question  of  required  sample  support,  a 
Monte  Garlo  analysis  was  performed.  This  analysis  showed  Reed’s  rule  does  hold  when  the  GSM 
algorithm  is  fully  adaptive.  However,  the  rule  does  not  hold  when  AToof  is  reduced  below  MN 
(fully  adaptive). 
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The  required  secondary  data  support  is  related  to  the  dimension  of  the  interference  subspace. 
The  Monte  Carlo  analysis  showed  iVooF  should  be  equal  to  the  dimension  of  the  interference 
subspace  for  optimal  SINK  performance.  For  the  interference  plus  noise  covariance  matrix  used 
in  this  experiment,  the  number  of  required  secondary  data  vectors  should  be  set  to  approximately 
2.5A^dof* 

As  a  final  point  of  interest,  the  GIF  non-homogeneity  detector  was  introduced  into  the  CSM 
algorithm.  A  worst  case  situation  was  used  in  that  the  non-homogeneity  inserted  into  the  data 
was  comprised  of  a  duplicate  target  in  another  range  cell.  This  represents  worst  case  because  the 
non-homogeneity  does  not  lie  in  the  interference  subspace.  The  CSM  algorithm  cannot  filter  it 
out.  The  GIF  non-homogeneity  detector  was  shown  to  significantly  improve  the  performance  of 
the  algorithm  in  this  case. 

The  most  significant  problem  with  implementing  the  GIF  non-homogeneity  detector  is  the 
requirement  of  a  fully  ranked  covariance  estimate  based  on  the  entire  secondary  data  set.  This 
requirement  conflicts  with  the  goal  of  reducing  the  necessary  secondary  data  support.  Using  the 
GIF  detection  scheme  with  the  CSM  algorithm  requires  a  minimum  of  MN  secondary  data  vectors. 
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VII.  Conclusions  and  Recommendations 


7. 1  Conclusions 

One  of  the  significant  design  considerations  for  any  STAP  algorithm  is  the  required  secondary 
data  support.  For  an  SMI  algorithm,  Reed  found  the  required  secondary  data  support  as  a  function 
of  the  number  of  DOF  [15]  used  in  the  adaptive  algorithm.  Through  an  exhaustive  derivation  based 
on  the  probability  of  detection  Pb)  Reed  determined  this  support  as  approximately  2iVboF‘  A 
Monte  Carlo  simulation  was  used  to  experimentally  determine  the  required  secondary  data  support 
for  the  CSM  algorithm. 

It  is  interesting  to  relate  Reed’s  rule  to  the  physical  size  of  the  interference  plus  noise  covari¬ 
ance  matrix  estimate  R.  For  SMI  techniques,  this  sample  matrix  is  immediately  inverted  to  obtain 
the  adaptive  weight  vector.  Hence,  it’s  dimension  is  N\)oy  x  ATdof-  The  number  of  secondary  data 
vectors  used  to  form  the  estimate  must  be  at  least  iVboF  or  the  matrix  is  singular  and  cannot  be 
inverted.  Using  2iVDOF  secondary  data  vectors  ensures  that  not  only  is  the  matrix  invertible  or 
non-singular  but  also  that  performance  of  the  STAP  algorithm  is  within  3  dB  of  optimal  SINR. 

The  CSM  algorithm  uses  the  eigenvalue/eigenvector  decomposition  of  a  full  dimensioned  R  to 
choose  ATdof  out  of  a  total  oiMN  DOF.  Correct  determination  of  the  eigenvalues  and  eigenvectors 
for  the  interference  and  noise  subspaces  require  a  full  rank  R.  However,  applying  Reed’s  rule  typ¬ 
ically  results  in  R  being  less  than  full  rank.  This  introduces  noise  into  the  eigenvalue/eigenvector 
decomposition.  This  calculation  noise  compromises  the  ability  of  the  algorithm  to  correctly  deter¬ 
mine  the  optimal  DOF  by  corrupting  the  metric  used  to  rank  order  the  DOF. 

The  problem  can  be  alleviated  through  a  careful  choice  of  iVooF-  Setting  this  value  equal 
to  the  dimension  of  the  interference  subspace  minimizes  the  impact  of  the  calculation  noise  in 
the  eigenvalue/eigenvector  decomposition.  The  minimization  occurs  because  the  calculation  errors 
in  the  eigenvalue/eigenvector  decomposition  have  greater  effect  on  the  DOF  corresponding  to  the 
noise  subspace  than  those  in  the  interference  subspace.  The  greater  effect  is  due  to  the  smaller 
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eigenvalues  associated  with  the  eigenvectors  spanning  the  noise  subspace.  Since  their  magnitude 
is  smaller,  the  relative  error  is  much  larger.  By  choosing  the  DOF  where  the  error  has  the  least 
impact,  the  efficiency  of  the  algorithm  is  optimized. 

Decreasing  A^dof  below  the  dimension  of  the  interference  subspace  results  in  a  corresponding 
decrease  in  SINK  performance.  This  decrease  in  performance  occurs  as  a  result  of  the  errors  in 
these  eigenvalues.  Although  the  relative  error  is  small  when  compared  to  the  magnitude  of  the 
eigenvalue,  it  is  enough  to  drastically  change  the  order  in  which  the  DOF  are  chosen.  Since  iVooF 
is  now  smaller  than  the  dimension  of  the  interference  subspace,  the  order  has  become  important. 
There  is  a  corresponding  decrease  in  performance. 

Increasing  iVboF  above  the  dimension  of  the  interference  subspace  also  results  in  a  decrease 
in  SINK  performance.  However,  the  algorithm  does  converge  to  optimal  as  secondary  data  set  size 
increases.  The  most  peculiar  behavior  occurs  when  the  secondary  data  set  size  approaches  MN^ 
the  full  rank  of  R.  There  is  a  large  drop  in  SINK  with  a  minimum  occurring  at  a  secondary  data 
set  size  J  of  MiV  —  1.  When  J  increases  beyond  MN  the  curve  resumes  the  expected  upward  slope 
towards  optimal  SINK. 

The  main  result  of  this  work  is  to  choose  iVboF  equal  to  the  dimension  of  the  interference 
subspace.  Obtaining  performance  within  3  dB  of  optimal  SINR  then  requires  2.5JVdof  secondary 
data  vectors.  Reed’s  rule  of  2iVDOF  was  found  to  only  apply  in  the  case  of  full  dimensionality.  The 
problem  with  this  approach  is  obtaining  an  accurate  estimate  of  the  interference  subspace  dimen¬ 
sion.  An  accurate  estimate  by  use  of  the  eigenvalues  or  singular  values  requires  more  secondary 
data  vectors  than  that  needed  to  reach  the  desired  SINR  performance.  A  method  of  determining 
interference  subspace  dimension  by  other  means  is  needed.  Brennan’s  rule  should  be  considered  for 
this  application. 

As  a  second  goal,  a  brief  analysis  of  non-homogeneities  in  the  secondary  data  was  performed. 
A  worse  case  scenario  for  the  CSM  algorithm  was  used  in  that  the  non-homogeneity  introduced  into 
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the  secondary  data  was  identical  to  the  desired  target  in  angle  and  Doppler.  The  non-homogeneity 
only  differed  in  range.  This  represents  worst  case  because  the  blocking  matrix  B  used  by  the  CSM 
algorithm  does  not  allow  this  signal  to  become  part  of  the  interference  plus  noise  estimate. 

The  CSM  algorithm  was  found  to  be  susceptible  to  the  single  non-homogeneity  in  the  sec¬ 
ondary  data.  The  difference  between  the  target  peak  and  the  next  highest  peak  in  the  algorithm’s 
output  was  reduced  by  approximately  3  dB.  Introducing  more  non-homogeneities  in  the  data  would 
typically  result  in  even  less  separation. 

Through  the  use  of  the  GIF  non-homogeneity  detector,  the  separation  between  target  and  next 
highest  peak  was  increased  by  3  dB.  This  effectively  negates  the  impact  of  the  non-homogeneity. 
The  GIF  works  very  well  with  the  CSM  algorithm  in  this  situation. 

7.5  Recommendations  for  Future  Research 

As  mentioned  in  Chapter  IV,  the  dimension  reducing  projection  matrix  U  in  the  CSM  al¬ 
gorithm  does  not  have  to  be  comprised  of  eigenvectors.  The  CSM  algorithm  provides  maximum 
SINK  only  for  the  case  when  eigenvectors  are  used  to  construct  this  projection  matrix. 

An  excellent  opportunity  for  future  research  is  development  of  a  different  projection  matrix. 
The  STAF  formulation  should  be  kept  in  the  GSC  architecture  because  it  offers  excellent  visual¬ 
ization  of  the  overall  process.  It  also  presents  a  good  analogy  to  the  Wiener  filter.  Avoiding  the 
eigenvalue/eigenvector  decomposition  would  result  in  several  advantages. 

If  a  constant  projection  matrix  can  be  defined,  the  calculation  costs  of  a  decomposition  can 
be  avoided.  This  cost  is  a  very  important  consideration  if  the  algorithm  were  to  be  considered  for 
real  time  implementation. 

A  second  advantage  to  a  constant  projection  matrix  is  secondary  data  support.  Avoiding  the 
decomposition  of  the  full  dimensioned  covariance  matrix  would  allow  the  algorithm  to  fall  under 
the  SMI  umbrella.  Hence,  Reed’s  rule  would  apply. 
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The  major  drawback  to  a  constant  projection  matrix  is  it  is  no  longer  adaptive.  The  same 
DOF  will  be  chosen  regardless  of  the  interference  environment  the  algorithm  is  placed  in.  Although 
the  algorithm  as  a  whole  will  remain  adaptive  due  to  the  reliance  on  the  interference  plus  noise 
covariance  estimate,  the  DOF  chosen  will  not  result  in  optimal  performance.  For  this  reason, 
a  constant  projection  matrix  may  not  be  feasible.  However,  it  is  important  to  consider  that  a 
constant  projection  matrix  forming  a  basis  for  likely  interference  scenarios  might  produce  acceptable 
performance,  albeit  sub-optimal  performance.  The  gains  in  calculation  cost  and  secondary  data 
support  from  avoiding  a  decomposition  might  far  outweigh  the  loss  in  performance. 

A  second  area  of  future  research  involves  examining  Brennan’s  rule  as  a  method  to  predict 
the  interference  subspace  dimension.  The  current  method  of  determining  this  dimension  involves 
the  use  of  either  the  eigenvalues  or  the  singular  values.  These  methods  typically  require  a  greater 
number  of  secondary  data  vectors  than  that  required  to  meet  the  desired  SINK  performance.  If 
Brennan’s  rule  gives  a  sufficiently  accurate  estimate,  the  algorithm  could  be  optimized  without  a 
penalty  in  secondary  data  support. 

Another  potential  area  for  research  is  the  GIF  non-homogeneity  detector.  As  mentioned 
previously,  the  major  drawback  to  this  method  is  the  inverse  of  a  covariance  matrix  estimate  R, 
based  on  the  full  secondary  data  set.  Performing  the  inverse  requires  the  estimate  be  full 
rank.  This  full  rank  requirement  then  translates  to  a  corresponding  requirement  of  secondary 
data  support.  In  the  case  of  the  GSM  algorithm,  choosing  Nbof  equal  to  the  dimension  of  the 
interference  subspace  means  only  2.5iVDOF  secondary  data  vectors  are  required.  However,  because 
R^  is  estimated  from  full  dimensioned  space- time  snapshots  there  are  at  least  MN  secondary  data 
vectors  required  to  ensure  it  is  invertible.  Typically,  2.5A^dof  <  MN,  A  new  method  of  detecting 
and  excising  non-homogeneities  would  alleviate  this  problem  of  secondary  data  support. 
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Appendix  A.  Kronecker  Product 


The  Kronecker  Product  is  simple  in  concept.  Given  two  matrices  A  and  B,  the  Kronecker  Product 
is  defined  in  the  following  manner  [9], 


A(giB 


.AiiB  i4i2B  A13B  •  •  •  AimB 

A21B  A22B  A23B  •  •  •  A2mB 

AjviB  Ajv2B  A7V3B  •  •  •  AjvmB 


(A.1) 


Therefore,  given  an  N  x  M  matrix  A  and  P  xQ  matrix  B  results  in  an  NP  x  MQ  matrix. 
The  following  table  of  identities  holds  for  Kronecker  products. 

Table  A.l  Properties  of  Kronecker  Products 

1.  (A  +  B)0C  =  A(g)C  +  B(8)C 

2.  (A  0  B)  0  C  =  A  0  (B  0  C) 

3.  0!(A  0  B)  =  (qjA)  0  B)  =  a  0  (ccB) 

4.  (A  0  B)^  =  A’’  0  B^ 

5.  (A  0  B)-i  =  A-i  0  B-i 

6.  (A0B)(C0D)  =  (AC)0(B0D) 

7.  (A  0  B)  =  (A  0 1)(I  0  B) 
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Appendix  B.  CSM  Blocking  Matrix  Determination 

Two  simple  methods  for  computing  the  signal  blocking  matrix  are  the  singular  value  decomposition 
algorithm  and  the  QR  decomposition  algorithm  [7].  The  singular  value  decomposition  of  a  matrix 
A  produces 

A  =  USV*^,  (B.l) 

where  the  matrices  U  and  V  are  unitary  and  S  is  diagonal.  The  left  and  right  singular  vectors 
are  contained  in  U  and  V,  respectively.  If  the  input  matrix  A  is  of  dimension  P  x  C,  then  U  is 
P  X  F,  V  is  C  X  C,  and  S  is  a  diagonal  P  xC  matrix  of  the  singular  values.  Note  that  regardless 
of  whether  the  input  matrix  is  real  or  complex,  the  nonzero  singular  values  are  always  real  and 
positive.  An  excellent  discussion  of  the  singular  value  decomposition  is  contained  in  [20], 

Using  these  definitions,  the  singular  value  decomposition  of  the  complex  conjugate  transposed 
steering  vector  (dimension  1  x  NM)  is 

=  USV^,  (B.2) 

where  U  is  1  x  1,  S  is  1  x  NM,  and  V  is  NM  x  NM.  Taking  the  hermitian  of  both  sides  results 
in 


s  =  VEU^.  (B.3) 

If  we  choose  the  blocking  matrix  such  that 

B  =  V(:,2:iVM)"  (B.4) 
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where  the  notation  used  is  Matlab®  notation  for  taking  all  the  rows  of  the  matrix  and  columns 
2  through  ATM,  then  the  product  of  the  blocking  matrix  and  steering  vector  is 

Bs  =  V(:,  2  :  (B.5) 


Let’s  examine  product  of  the  reduced  and  full  dimension  right  singular  vector  matrices.  First 
express  them  in  terms  of  the  right  singular  vectors  themselves, 

Vl  V2  •  •  •  vatm 

where  the  last  step  takes  advantage  of  the  orthonormal  columns  of  V. 

Since  the  input  matrix  to  the  singular  value  decomposition  is  actually  a  row  vector,  there  will 
be  only  one  singular  value  for  the  decomposition  of  s^,  i.e. 


0  0 

0  Inm-1 


(B.6) 


V(:,2  :iV’M)^V  = 


vf 


V3^ 


^NM 


s  = 


(Tl  0  •  •  • 


(B.7) 


Recognizing  this  fact  and  substituting  Eqn.  (B.6)  into  Eqn.  (B.5)  results  in 


Bs  = 


0  0 
0 


SU"  =  0, 


the  desired  result. 


(B.8) 
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The  second  method  mentioned  by  Goldstein  is  the  QR  decomposition.  This  decomposition 
produces 

A  =  QR,  (B.9) 

where  Q  is  P  x  P  and  R  is  an  upper  triangular  matrix  of  dimension  P  xC  given  A  is  P  x  C7.  The 
elements  of  R  are  given  by  [20] 

Tij  =  qf  a,-,  (B.IO) 

where  q*  and  a.j  are  the  columns  of  Q  and  A,  respectively. 

The  matrix  Q  is  calculated  using  the  Gram-Schmidt  orthogonalization  process,  hence  the 
columns  of  Q  are  orthonormal  and  the  matrix  is  termed  orthogonal.  The  orthogonalization  process 
produces  a  set  of  vectors  that  when  taken  as  a  whole  represent  the  subspace  spanned  by  the  input 
matrix. 

The  QR  decomposition  of  the  steering  vector  is 

s  =  QR,  (B.ll) 

where  Q  is  NM  x  NM  and  R  is  iVM  x  1.  The  only  nonzero  element  contained  in  R  is  the  first 
element.  Therefore,  choosing  the  blocking  matrix  as 

B  =  Q(;,2:i\rM)"  (B.12) 

results  in  Bs  =  0  as  desired.  This  can  be  shown  in  the  exact  same  manner  as  that  done  for  the 
singular  value  decomposition. 
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Appendix  C,  Performance  Measures 


The  main  goal  of  the  radar  is  to  detect  the  target  amid  noise,  interference,  and  jamming.  Target 
detection  is  a  binary  decision,  either  a  target  is  present  or  there  is  not.  Therefore,  one  of  the 
obvious  metrics  for  algorithm  comparison  is  the  decision  test.  This  section  addresses  several  of  the 
common  test  statistics  used  for  the  binary  detection  problem,  beginning  with  the  one  developed  by 
Brennan  and  Reed  in  their  pivotal  paper  [2].  Then  the  discussion  transitions  to  two  test  statistics 
offering  the  desirable  constant  false  alarm  rate  (CFAR)  characteristic. 

CA  RMB  Test 

Continuing  exactly  as  in  Brennan  and  Reed’s  earlier  work  [2],  the  radar  is  cast  in  the  light  of 
a  binary  hypothesis  problem.  Either  there  is  a  target  in  the  data  and  the  hypothesis  is  Hi  or  there 
is  not  and  the  hypothesis  is  The  data  models  are  expressed  by 


ffo: 

Xt=  ^ 

Hi  : 

Xi=  as  +  n, 

(C.l) 

where  n  is  the  colored  noise  with  associated  covariance  matrix  estimate  R/,  s  is  the  normalized 
steering  vector,  and  a  is  a  random  complex  coefficient  whose  phase  (j>  is  uniformly  distributed 
between  0  and  27r.  The  subscripted  I  emphasizes  the  fact  that  the  data  vector  under  consideration 
corresponds  to  a  particular  range  cell.  This  is  the  exact  same  method  used  by  Reed  to  develop 
the  adaptive  processor  weights  that  maximized  the  probability  of  detection  and  subsequently  the 
SNR  [2]  (introduced  in  Section  4.2). 
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Under  these  hypotheses  and  given  complex  Gaussian  noise,  the  probability  density  functions 
of  the  incoming  data  vector  Xi  conditioned  on  the  random  phase  of  the  target  steering  vector  are 


fxi\Hoy4> 


_ I _ _ _ ^-(Xi-olb)^Ai  ^(Xi-as) 

^JVMdet(R/) 

_ L^e-xr^rx, 

^ArMdet(R;)  ’ 


(C.2) 


where  N  is  the  number  of  spatial  degrees  of  freedom  and  M  is  the  number  of  temporal  degrees  of 
freedom  such  that  the  steering  vectors  a  and  b  are  iV  x  1  and  M  x  1  respectively.  This  produces 
a  space-time  steering  vector,  v  or  s,  of  size  NM  x  1.  The  likelihood  ratio  test  then  becomes, 


Ai 


Hi 

^0  Vu 


(C.3) 


where  £  is  the  expected  value  operator  with  respect  to  the  density  of  the  random  phase  and 
7]  is  the  threshold.  Using  the  densities  given  above,  the  likelihood  ratio  test  is 

Hi 

=  (C.4) 


where  loix)  is  the  modified  Bessel  function  of  the  first  kind.  This  is  a  Bayes  optimal  test.  Isolating 
the  signal  power  expression  on  the  left  side  results  in  the  Reed,  Mallett,  and  Brennan  (RMB) 
test  [2], 


Hi 


Ai  =  |s"R,  xif 


2|a| 


=  771. 


(C.5) 


Since  the  new  threshold  is  a  function  of  the  covariance  matrix  and  incoming  data,  it  cannot  be 
determined  beforehand  to  achieve  a  specified  false  alarm  probability.  Hence,  the  test  is  not  CFAR 
and  is  not  practical  for  use  as  a  test  statistic. 
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a2  GLRTTest 


The  Generalized  Likelihood  Ratio  Test  (GLRT)  was  developed  by  Kelly  in  a  paper  on  adaptive 
detection  [13].  The  test  is  given  by 


(s^R-r^s)  [i + Rr^xj)] 


Hi 

^0 


(C.6) 


where  R/  represents  the  interference  covariance  matrix  estimate,  Xi  Ihe  primary  data  vector,  s  the 
steering  vector,  and  2J  is  the  number  of  secondary  data  vectors  used  to  calculate  the  covariance 
matrix  estimate.  Kelly  shows  this  test  does  in  fact  possess  the  generalized  CFAR  property  via  the 
normalization  by  the  signal  dependent  quantity  in  the  denominator. 


as  AMF  CFAR  Test 


The  final  CFAR  test  under  discussion  is  one  introduced  just  recently  (1991)  by  Chen.  The 
adaptive  matched  filter  CFAR  (AMF  CFAR)  test  [4]  is 


As  = 


|s"Rr^x<P 

s^Rr^s 


Hi 

^0 


(C.7) 


This  test  statistic  has  also  been  referred  to  as  the  Nusui  test  statistic  in  articles  by  Melvin.  Because 
the  test  provides  a  constant  threshold  that  is  a  function  of  the  probability  of  false  alarm,  it  has  an 
embedded  CFAR  characteristic.  This  is  the  source  of  the  name. 

The  relationship  between  the  threshold  and  PpA  is  given  by  [4] 


Pfa  = 


(C.8) 


where  773  is  the  threshold. 

This  test  was  developed  by  Chen  and  Reed  as  an  improvement  of  the  GLRT  test. 
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model,  4-5 
output  (j/),  4-12 
weight  vector  (wb),  4-10 
weight  vector  (wp),  4-12 
weight  vector  (wqsc),  4-8 
operation  theory,  4-10 
reduced  dimension,  4-4,  4-11 
block  diagram,  4-14 
output  (2/r),  4-13 
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projection  matrix,  4-13 
weight  vector  (w*),  4-15 
transformation  matrix  (T),  4-6 
GLRT  test,  see  threshold  detection,  GLRT 
test 

grating  lobes,  2-17 
grazing  angle  3-23 

Hadamard  product  (©),  3-41 
Hamming  window 

specifications,  3-55 
Hanning  window 

specifications,  6-9 

interference  covariance,  see  covariance 
matrix  (R) 

Jammer-to-Noise  Ratio  (JNR),  3-17 

Kronecker  delta  (5),  3-14 
Kronecker  product  (©),  3-11,  A-1 

likelihood  ratio  test,  C-2 

Maximum  Likelihood  Estimation  (MLE),  see 
covariance  matrix  (R),  estimation 
sample  support  for,  6-1 
MCARM,  1-3,  5-7 
Monte  Carlo  simulation,  6-1,  6-11 


noise  covariance,  see  covariance  matrix  (R) 
noise  model,  see  STAR,  model,  thermal  noise 
noise  power  per  element  (cr^),  3-14 
non-homogeneities 
definition  of,  5-5 
detection  of,  5-4 

plane  wave,  2-1 

frequency  spectrum  of,  2-3 
principal  components,  4-29 
method  of,  4-15 
probability  of  detection,  4-2 
maximized,  4-2 

probability  of  false  alarm,  4-1,  5-13,  C-3 
propagation  velocity  (t;p),  2-5,  2-8,  2-10,  3-7 
Pulse  Repetition  Interval  (PRI),  3-4 

Q-function,  4-2 
QR  decomposition,  B-3 

Radar  Cross  Section  (RCS) 
of  clutter  patch,  3-25 
receiver  losses  (Lr),  3-18 
Reed’s  rule,  5-3 

in  CSM  algorithm,  6-1,  6-9,  6-11,  6-13, 
6-17,  6-22 

RMB  test,  see  threshold  detection,  RMB  test 
Sample  Matrix  Inversion  (SMI),  5-3,  6-1 
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sampling 


space-time  signal,  2-4,  2-13,  2-14 


of  clutter,  3-31 
spatial,  2-6 
temporal,  2-7,  3-9 
secondary  data,  5-2,  5-13 

analysis  of  required  support,  6-9,  6-22 
Singular  Value  Decomposition  (SVD),  6-19, 
B-1 

SINK,  3-41 

for  CSM  algorithm 
fully  adaptive,  6-2 
partially  adaptive,  6-6 
for  FTS,  3-53 

for  tapered  fully  adaptive  STAP,  3-42 
maximized,  3-40 

optimum  (fully  adaptive  STAP),  3-41 
using  interference  covariance  matrix  esti¬ 
mate  (R),  5-3 

SINK  Loss  (Lsinr),  3-43,  3-53 
slowness  vector  (a),  2-2 
SNR 

loss  due  to  interference  covariance  matrix 
estimation,  5-3 
obtaining  maximum,  4-3,  4-5 
optimum,  3-44 

Space-Time  Adaptive  Processing,  see  STAP 


superposition  of  plane  waves,  2-4 
space-time  snapshot  (x),  3-12 
in  clutter  model,  3-24 
in  jammer  model,  3-17 
in  noise  model,  3-14 

space-time  steering  vector,  see  steering  vector 
spatial  frequency,  see  frequency,  spatial 
STAP,  3-1 

block  diagram,  3-44 
CSM  model,  see  Generalized  Sidelobe 
Canceler  (GSC) 

FTS,  see  factored  time-space  (FTS) 
fully  adaptive,  3-38 
weight  vector  (w),  3-40 
GSC  form,  see  Generalized  Sidelobe 
Canceler  (GSC) 
model 

clutter,  3-21 

electromagnetic  environment,  3-2,  3-4 
jamming,  3-14 
physical  environment,  3-2 
thermal  noise,  3-13 
steering  vector 

in  jammer  model,  3-16,  3-20 
normalized  space-time  (s),  4-4 
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space-time  (v),  3-11 
spatial  (a),  3-10 
temporal  (b),  3-11 

temporal  frequency,  see  frequency,  temporal 
threshold  detection 

AMF  CFAR  test,  C-3 
GLRT  test,  C-3 
RMB  test,  C-2 

wavelength  (A),  2-2,  2-10 
wavenumber  vector  (k),  2-3 
wavenumber-frequency  space,  2-3,  2-4 
weight  vector  (w) 

describing  entire  CSM  algorithm 
fully  adaptive,  6-3 
partially  adaptive,  6-7 
tapered,  3-40 
whitening  filter,  5-6 
windowing,  2-17 

Blackman  Harris,  see  Blackman  Harris 
window 

Hamming,  see  Hamming  window 
Hanning,  see  Hanning  window 
in  Doppler  filter,  3-51 
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